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ABSTRACT 


The  determination  of  the  Mode  I  stress  intensity  factors 
for  selected  crack  configurations,  using  finite  element  methods  and 
energy  release  rate  principles,  is  the  subject  of  this  study.  The 
crack  configurations  which  were  investigated  were  the  double  edge 
crack,  the  single  edge  crack,  and  the  center  crack.  The  method  of 
analysis  utilized  was  the  "Stiffness  Derivative  Method”  [7].  This 
approach  relates  the  change  in  strain  energy  resulting  from  crack 
advancement,  to  the  change  in  the  stiffness  matrix  of  the  struc 
ture  containing  the  crack.  The  results  indicated  that  through 
mesh  optimization  and  proper  control  of  certain  parameters 
including  the  crack  advance  increment,  the  crack  tip  element  con 
tour  size,  and  mesh  refinement,  an  accurate  solution  can  be  calcu-- 
lated  with  a  relatively  coarse  finite  element  mesh  consisting 
entirely  of  contemporary  elements.  The  numerically  generated 
solutions  are  compared  with  analytical  solutions  with  the  results 
within  0.001  percent  of  each  other  for  the  double  edge  crack, 

0.858  percent  for  the  single  edge  crack,  and  2.021  percent  for  the 


center  crack. 
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NOMENCLATURE 


Symbol 


Meaning 


[B] 

b 

[D] 

E 

e 

F 

G 

^1,  KlI, 
Kin 

[K] 

[k] 

{L} 

O 

P 

r 

S 


u,v,w 


Crack  half-length 

Strain-displacement  proportionality  matrix 
Plate  half-width 

Stress-strain  proportionality  matrix 
Young’s  modulus 
Primary  subscript  designator 
Work  done  by  external  forces 
Energy  release  rate 

Stress  intensity  factor  (Modes  I,  II,  III) 

Overall  stiffness  matrix 
Element  stiffness  matrix 
Load  vector 

Element  superscript  designator 

Secondary  subscript  designator 

Total  potential  energy 

Radial  distance  from  crack  tip 

Thermal  and  residual  effect  energy  term 

X,  y  and  z  displacements 


Vlll 


{u} 

U 

W 


T 


0 
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Nodal  displacement  vector 

Elastic  strain  energy 

Energy  available  for  crack  growth 

Normal  stress  component 

Stress  vector 

Shear  stress  component 

Strain  vector 

Poisson^ s  ratio 

Angle  of  inclination  from  crack  line 
Dummy  integration  variable 
Crack  tip  closure  increment 
Strain  energy  density 
Crack  tip  inner  element  contour 
Crack  tip  outer  element  contour 


INTRODUCTION 


The  initiation  and  growth  of  cracks  in  a  material  as  a 
result  of  pre-existing  flaws  has  been  an  extremely  important  phe¬ 
nomenon  in  the  historical  development  and  use  of  materials.  With 
the  advent  of  today's  high-strength  materials  and  their  subsequent 

application  in  many  sectors  of  our  society,  a  more  definitive 
understanding  of  the  fracture  phenomenon  is  essential  if  these 

materials  are  to  be  utilized  safely  and  most  effectively.  This 
need  is  readily  apparent  in  many  of  the  current  high-technology 
design  applications  which  employ  advanced  composites  with  spe 
cified  material  properties  and  a  high  strength-to-weight  ratio. 

As  a  result  of  such  advanced  application  it  is  no  longer  practical 
in  design  to  overspecify  the  material  with  a  sufficiently  large 
factor  of  safety  to  compensate  for  any  probable  defects.  Instead, 
the  designer  must  work  within  increasingly  smaller  safety  factors 
while  accounting  for  the  inherent  defects  present  in  the  material 
resulting  from  the  fabrication  process  as  well  as  from  damage 
incurred  subsequent  to  fabrication. 

To  define  the  stress  field  around  a  specific  crack  con¬ 
figuration,  it  is  necessary  to  determine  the  stress  intensity  fac¬ 
tor.  With  that  determined,  one  is  then  able  to  predict  the  crack 
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growth  and  the  time  to  fracture  which  are  the  important  factors 
related  to  any  defect  [1].  Of  the  three  crack  mode  stress  inten¬ 
sity  factors,  the  Mode  I  factor,  Kj,  is  generally  the  most  impor¬ 
tant,  and  is  the  one  which  is  evaluated  in  this  report  [!]• 

The  modeling  of  the  crack  configurations  and  determination 
of  the  respective  stress  intensity  factors  is  accomplished  by 
means  of  numerical  methods;  specifically  the  finite  element 
method.  Some  previous  studies  using  finite  elements  to  solve  this 
problem  have  developed  special  "singularity  elements”  for  use  at 
the  crack  tip  [2-4]  .  Such  elements  are  specifically  designed  to 
incorporate  the  elastic  stress  singularity  at  the  crack  tip.  The 
method  is  classified  as  one  of  the  direct  methods.  While  reaso¬ 
nable  accuracy  was  obtained  by  these  means,  it  is  more  desirable 
to  use  conventional  elements  in  view  of  their  universal  availabi¬ 
lity  and  less  complicated  incorporation  into  the  overall  mesh. 

In  selecting  a  method  of  analysis  as  a  basis  for  this 
study,  a  number  of  factors  serve  as  guidelines.  First,  the  model 
must  use  conventional  elements  for  reasons  previously  discussed. 

In  addition,  error  based  on  comparisons  with  analytical  solutions 
should  initially  be  within  5  percent  to  warrant  further  investiga¬ 
tion  and  refinement.  The  term  "error",  as  used  in  this  report,  is 


defined  as: 


3 


(Analytical  solution)^(Numerical  solution) ^ 
(Analytical  solution)  _j 

The  third  guideline  is  that  the  mesh  must  be  optimized  so  as  to 
result  in  a  relatively  coarse  mesh  configuration  which  is  still 
sufficiently  accurate.  The  mesh  optimization  is  important  so  as 
to  minimize  computing  time  and  thereby  lay  the  groundwork  for  the 
implementation  of  the  method  in  solving  the  more  extensive  three- 
dimensional  problems  related  to  interlaminar  flaws  in  composite 
laminates.  Based  on  the  specified  criteria,  a  number  of  the 
indirect  methods  involving  energy  methods  have  been  found  to  be 
viable  choices  [5-7].  Of  those,  the  "Stiffness  Derivative  Method" 
used  by  Parks  [7]  has  been  chosen  for  this  study. 

The  "Stiffness  Derivative  Method”  is  based  on  evaluating 
the  energy  release  due  to  crack  extension  as  it  relates  to  the 
corresponding  change  in  the  stiffness.  It  requires  the  results  of 
only  one  crack  advancement  to  compute  the  solution. 

In  this  current  research  endeavor,  the  method  is  used  to 
evaluate  the  Mode  I  stress  intensity  factors  of  three  two- 
dimensional  crack  configurations.  Those  configurations  include 
the  single  edge  crack,  the  double  edge  crack,  and  the  center 
crack.  In  all  cases,  mesh  optimization  techniques  [8]  are 
employed,  and  numerical  solutions  generated  which  are  compared  to 
their  corresponding  analytical  solutions. 


Percent  Error 


CHAPTER  1;  THEORETICAL  FORMULATION 


1 •!  The  Stress  Intensity  Factor 

The  characterization  of  the  deformation  of  a  crack  can  be 
separated  into  three  distinct  modes  of  local  deformation  as  shown 
in  Figure  1.  In  each  case,  the  crack  deformation  is  the  result  of 
a  different  state  of  stress  being  present  in  the  body.  The 
opening  mode  (Mode  I)  is  the  result  of  an  in-plane  normal  stress 
on  the  crack.  The  sliding  mode  (Mode  II)  also  involves  an  in¬ 
plane  stress,  but  in  this  case  it  is  a  shear  stress.  By  contrast 
the  presence  of  an  out-of-plane  shear  stress  causes  the  tearing 
mode  (Mode  III)  .  Very  often,  the  stress  state  near  the  crack  tip 
is  actually  a  mixed  mode  situation  with  a  combined  effect  of  the 
three  modes  to  varying  degrees.  Of  the  three  modes.  Mode  I  has 
been  found  to  be  the  predominant  deformation  mode,  and  is  the  one 
which  is  concentrated  on  in  this  study  [1]. 

In  the  case  of  an  infinite  plate  in  the  x-y  plane  with  a 
through  crack  and  an  in-plane  stress  state  normal  to  the  crack 
(Mode  I),  the  components  of  the  elastic  stress  field  of  an  element 
dxdy  at  a  distance  r  from  the  crack  tip  and  an  angle  0  from  the 
crack  line,  as  shown  in  Figure  2,  are  found  to  be  (omitting  higher 
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Figure  1.  Crack  Deformation  Modes 
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r 


Figure  2 


Crack  Tip  Stress  Field  Components  and  Coordinates 
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order  terms)  [1]: 


ax  =  a(a/2r)^2  cos(e/2)[l  -  sin(0/2)sin(3e/2) ]  (1.1) 

Oy  =  a(a/2r)^2  cos(0/2)[l  +  sin(0/2)sin(30/2)]  (1.2) 

Txy  =  a(a/2r)^2  sin(0/2)cos(0/2)cos(30/2)  (1.3) 

©2=0  (plane  stress)  (1.^) 


~  ^y)  (plane  strain) 


(1-5) 


In  these  expressions,  a  is  the  half~crack  length*  These  equations 
can  generally  be  expressed  as  [1]: 


(2TTr)V2 


where 


Ki  =  a(Tra)^'^ 


(1.6) 


(1.7) 


Kl  is  called  the  Mode  I  stress  intensity  factor. 

The  Mode  I  elastic  stress  field  of  the  element  dxdy  and 
the  Cartesian  displacements  (u,v,w)  can  be  expressed  as  [9]: 

Ox  =  cos(0/2)[l  -  sin(0/2)sin(30/2)] 

(2Trr)V2 


(1.8) 
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For  plane  strain: 
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~  ®y)  (1*19) 

Txy  =  Tyz  =  0  (1.20) 

u  =  2Kii(1+v)  (r/2Tr)V2  sin(9/2)[2  -  2v  +  cos^(0/2)]  (1.21) 
E 

V  =  i^^Ilill^(r/2x)l/2  cos(6/2)[-l  +  2v  +  sin^(9/2)]  (1.22) 
E 

w  =  0  (1.23) 

The  Mode  III  elastic  stress  field  and  displacements  are 
expressed  as  [9]: 


txz  =  — ^ sin(9/2) 

(2iTr)l/2 

(1.24) 

Ty2  =  cos(0/2) 

(2iTr)V2 

(1.25) 

'Jx=  Oy  =  Oz  -  Txy  ~  0 

(1.26) 

„  «  2Kiii(1+v)  ^2r/Tr')l/2  sinf9/2) 

(1.27) 

E 

o 

II 

> 

II 

(1.28) 

In  these  equations,  Kji  and  Kjn  are  the  Mode  II  and  Mode  III 
stress  intensity  factors  respectively. 

The  expression  for  Kj  in  equation  (1*7)  is  that  corres¬ 
ponding  to  the  case  of  a  crack  in  an  infinite  plate.  In  actual 
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applications  with  plates  of  finite  width,  the  expression  becomes 
less  accurate  as  the  ratio  of  the  crack  half-length  to  plate  half¬ 
width,  (a/b),  increases.  To  account  for  such  effects,  the  stress 
intensity  factor  is  expressed  as  [1]: 

Kj  =  a(TTa)^2  f(a/b)  (1.29) 

in  which  the  function  f(a/b)  is  a  function  of  the  crack  half- 
length,  a,  and  the  plate  half-width,  b. 

If  the  applied  stress,  a,  on  the  plate  with  the  crack  is 
increased,  it  will  eventually  reach  a  point  at  which  time  the 
plate  will  fracture.  The  applied  stress  at  this  point  is  the 
failure  stress,  Oq,  The  corresponding  stress  intensity  factor  is 
the  critical  stress  intensity  factor  or  the  fracture  toughness 
which  is  expressed  as  [1]: 

Ki  =:  (ira)^2  f(a/b)  (1.30) 

c 

Physically,  the  fracture  toughness  is  an  indication  of  the 
resistance  of  a  material  containing  a  crack  to  fracture.  With 
this  factor  known  for  a  certain  material,  it  is  possible  to  deter¬ 
mine  the  maximum  allowable  crack  size  without  the  onset  of 
fracture  for  a  given  stress  level.  Conversely,  the  fracture 
strength  of  the  material  for  a  given  crack  size  can  also  be  calcu¬ 
lated.  The  results  of  equation  (1.30)  are  valid  when  a  condition 
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of  plane  strain  exists.  If  a  plane  stress  condition  is  present, 
then  the  critical  stress  intensity  factor  is  dependent  upon  the 
thickness  of  the  plate  also  [1].  The  importance  of  knowing  the 
stress  intensity  factor  is  quite  apparent  since  with  that  infor¬ 
mation,  one  can  determine  the  crack  growth  behavior  as  well  as  the 
point  at  which  a  material  will  fracture. 


1 .2  The  Griffith  Energy  Criterion 

For  an  infinite  plate  with  a  through  crack  of  half-length, 
a,  and  normal  in-plane  stress,  o,  Griffith  showed  that  in  order 
for  the  crack  to  propagate,  the  following  conditions  must  exist 
[1]: 


(1.31) 

In  this  equation,  W  is  the  energy  available  for  crack  growth,  F  is 
the  work  done  by  external  forces,  and  U  is  elastic  strain  energy 
of  the  plate.  In  essence,  the  equation  indicates  that  crack  pro¬ 
pagation  will  occur  when  conditions  are  such  that  the  energy 
released  as  a  result  of  the  crack  extension  is  equal  to  the  energy 
required  for  the  crack  growth. 


If  the  loading  of  the  plate  is  such  that  the  ends  are 
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fixed,  then  the  external  forces  do  no  work.  In  this  situation, 
the  energy  for  crack  growth  comes  from  a  release  of  elastic 
energy,  resulting  in  a  decrease  in  the  elastic  energy  content  of 
the  plate.  If  the  ends  of  the  plate  are  free,  then  the  external 
forces  do  work  in  this  case.  The  result  of  the  input  work  energy 
is  an  increase  in  the  elastic  energy  content  of  the  plate,  as  well 
as  energy  being  made  available  for  crack  growth. 


The  concept  of  the  energy  release  rate,  as  developed  by 
Irwin  [10],  can  be  used  to  relate  the  stress  intensity  factor  to 
the  Griffith  criterion  for  crack  growth.  Physically,  the  energy 
release  rate,  G,  is  the  rate  at  which  energy  is  made  available  for 
crack  propagation,  with  the  energy  originating  from  a  release  of 
elastic  energy  of  the  material,  or  from  work  done  on  the  body  by 
external  forces,  or  some  combination  of  both.  The  energy  release 
rate  is  expressed  by  Irwin  in  terms  of  the  crack  closure  integral 
[9]; 


G 


lim 

6^0 


dx 


(1.32) 


This  is  an  expression  of  the  rate  of  change  of  work  energy  with 
respect  to  crack  length  as  a  result  of  closing  the  crack  tip  a 
distance  <S ,  as  in  Figure  3.  The  stresses  and  displacements  are 
those  of  the  crack  surface  which  is  pulled  closed.  In  equation 
(1.32),  the  factor  of  2  is  present  because  two  crack  faces  are 
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Figure  3.  Crack  Tip  Closure 
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being  moved  in  the  closure  process*  The  divisor  of  2  in  each  term 
in  the  integral  results  from  the  fact  that  the  stresses  which 
close  the  crack  are  zero  at  the  tip  and  increase  moving  further 
along  the  crack.  Hence,  the  2  in  the  denominator  is  related  to 
the  mean  stress  applied  to  close  the  crack.  The  stresses  are  eva¬ 
luated  at  r  =  X  and  9=0.  The  displacements  are  calculated  using 
r  =  5-x  and  0  =  tt  . 


Evaluating  the  first  term  of  the  integral,  the  stress  and 
displacement  terms  are  those  of  the  Mode  I  state  as  found  in 
equations  (1.9)  and  (1.14)  respectively.  Substituting  the  values 

for  r  and  0 : 

(1-33) 

(2ttx)1/2 


V 


4Ki(1-v^) 

E 


(1.34) 


The  integral  then  becomes: 


lira  2Ki^(I-v^) 
airE 


(1.35) 


Gi  is  the  Mode  I  component  of  the  total  energy  release  rate  which 
is : 


G  =  Gi  +  Gii  +  Gin 


(1.36) 


The  substitution  of  x/6  =  sin^(|>  in  equation  (1.35)  results  in: 
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2  2  .  2^  -,Vl 

_  2Kx  (1-v  )  r  r - iifL-SL  I  (26  )sin(l)cos<{i  d<t> 

'5^0  5  L  sln^il)  -* 

(1.37) 

=  ’'/^cos%  d<!)  (1.38) 

6+0  j.  ^  ^ 

ttE  0 


2  2  .  1  — /  2 

lim  4Kx  (l-v  )  r <1^  .  sin2j)  I  n,39) 

s-o  — 12  *  Jo 


=  )  (Plane  Strain)  (1.40) 

E 


For  plane  stress: 


Gx  =  Ki^/E 


(1.41) 


The  evaluation  of  the  second  terra  in  the  integral  expression  of 
equation  (1.32)  is  made  using  the  stress  and  displacement  ex¬ 
pressions  of  equation  (1.18)  and  (1.21)  respectively.  The  cal¬ 
culations  are  the  same  as  those  of  Gx  except  for  the  presence  of 
Kii  in  the  equations  instead  of  Kx*  The  result  is  the  Mode  II 
component  of  the  energy  release  rate: 

Gxx  =  (Plane  Strain)  (1.42) 

E 


or: 


(Plane  Stress)  (1*43) 


The  third  terra  in  the  integral  expression  of  equation 


(1.32)  is  calculated  using  equations  (1.25)  and  (1.27)  for  the 
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stress  and  displacement  terms  respectively.  The  procedure  is  the 
same,  with  a  difference  arising  in  the  resultant  form  of  the 
displacement  expression  as  compared  to  the  previous  calculations. 
In  this  case; 


= 

(2Trx)l/2 


(1.44) 


=  (2r/ir)l/!2 

E 


(1.45) 


The  Mode  III  component  of  the  energy  release  rate  is: 


lira  2Kiii  (1+v) 
— 


6-x  1 1/2 


/  [  —  ] 
J)  L  X  J 


(1.46) 


which  reduces  to: 

Gin  =  ^ni^itl)_  (1.47) 

E 

The  evaluation  of  equation  (1#32)  is  at  this  point 
complete,  with  G  representing  the  total  energy  release  rate.  The 
Mode  I  component,  for  the  plane  strain  case  is  the  component 

that  is  directly  applicable  to  this  research  effort,  and  which  is 
utilized  in  subsequent  computations. 
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1.3  The  Stiffness  Derivative  Development 


The  stiffness  derivative  method  as  outlined  by  Parks  [7] 
is  based  upon  the  relationship  between  the  change  in  total  poten¬ 
tial  energy  of  a  linear  elastic  body  with  a  crack,  and  the  change 

in  the  element  stiffness  matrices  of  the  finite  element  solution 
of  the  problem,  both  of  which  result  from  advancement  of  the 

crack.  The  total  potential  energy  expression  as  it  relates  to  the 
finite  element  solution  and  the  stress  intensity  factor  is  devel¬ 
oped  in  this  work  based  on  work  done  by  Hellen  [11]  and  Parks 
[7].  Development  of  the  elasticity  equations  for  the  finite 
element  solutions  involves  the  use  of  matrix  notation  for  the 
variables,  as  the  equations  are  representative  of  systems  with 
many  degrees  of  freedom.  For  example,  the  finite  element  solutions 
of  this  work  employ  the  two-dimensional,  four-node  quadrilateral 
element.  The  element  has  two  degrees  of  freedom  at  each  node  and 
its  nodal  displacement  vector  is  expressed  as: 


{u}^  =  (uj^,v2^,U2,V2,U2,V2,u^,v^) 


(1.48) 


The  I  denotes  that  this  is  the  vector  for  an  individual  element, 
and  T  specifies  the  transpose  of  the  designated  matrix.  The  total 
strain  vector  at  a  point  on  the  element  is: 


{£}  =  {e}e  +  {£}. 


(1.49) 
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in  which  {e}g  represents  the  primary  strains  which  result  from 
mechanical  forces,  and  {e}©  represents  the  secondary  strains  which 
are  due  to  thermal  or  residual  effects.  The  stress  vector  at  this 
point  is: 

{a}  =  [D]{e}e  =  [D]({e}  -  {e}o)  (1-50) 

In  this  equation,  [D]  is  a  proportionality  matrix  consisting  of 
elastic  constants  of  Young's  modulus  E  and  Poisson's  ratio  v.  It 
is  a  square  symmetric  matrix,  and  therefore: 

[D]T  =  [D]  (1«51) 

The  total  strain  vector  and  the  nodal  displacement  vector  are 
related  by  means  of  the  [B]  matrix,  the  terms  of  which  are  deriva¬ 
tives  of  the  element  shape  function  at  the  specified  point: 

{e}  =  [B]{u}Ji  (1.52) 

The  expression  for  the  strain  energy  density  at  the  point  of 
interest  is: 

£2  =l/2{a}T{e}e  =  l^{a}T({e}  -  {e}^)  (1.53) 

From  equations  (1.50)  and  (1.51): 

{a}T  =  ({e}T  -  {£}o)[D]T  =  ({£}T  -  {e}Q)[D]  (1.54) 


From  equation  (1.52): 
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{e}^  =  {u}^T[b]T  (1.55) 

Substituting  the  previous  two  equations  into  equation  (1.53) 
results  in: 

n  =  l/2({u}^^[B]^  -  {e}o)  [D]  ([B]{u}^  -  {e}^)  (1-56) 

which  expands  to: 

=  V^{u}JIT[B]T[D][B]{u}A  -  iy^{u}AT[B]T[D]{e}Q 

-  l/2{e}T[D][B]{u}^  +  ]y2{£}T[D]  (1.57) 

The  third  term  of  this  expression  is  equal  to  its  transpose,  and 
can  therefore  be  combined  with  the  second  term  to  produce: 

n  =  V2{u}^T[b]T[d][B]{u}^  -  {u}^T[b]T[D]{£}^ 

+  V2{e}T[D]{£}^  (1-58) 

The  strain  energy  of  the  element  is  calculated  by  integrating 
the  strain  energy  density  of  the  element  volume  V: 

=  /  MV  =  1/2{u}^T  [/[B]T[D][B]dv]{u}£ 

I  I 

-  {u}AT  /[3]T[D]{£}  dV  +  V2  /{£}T[D]{£}^dV  (1.59) 

I  A  ° 


The  total  potential  energy  of  the  element  is  defined  as: 
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p£  =  _  {u}^'^{L}e 


(1.60) 


where  {L}e  Is  the  mechanical  load  vector  which  causes  the  primary 


strains.  In  addition; 


=  /[B]T[D]{e}  dV 
I 


(1.61) 


where  {L}©  is  the  load  vector  related  to  the  secondary  strains  in 
the  element.  The  element  stiffness  matrix  is  defined  as: 


[K]A  =  /[B]T[D][B]dV 
I 


(1.62) 


=  l^/{e}T[D]{e}QdV 


(1.63) 


Using  equations  (1.61)  through  (1.63),  the  total  potential  energy 
can  be  rewritten  as: 


n  -  l/5{u|‘^[Kl(„}^  -  |»}«(L|‘  .  s*  -  (u)‘T|l.)‘  (1.64) 

This  expression  is  defined  for  the  overall  structural  application 
by  deleting  the  superscripts  Z  designating  element  application. 

The  total  load  vector  can  be  written  as: 


ZT!r\<^ 


(1.65) 


In  the  overall  sense,  the  load  vector  {L}e  corresponds  to  loads 
that  are  applied  externally  since  the  internal  forces  all  cancel 
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at  the  internal  nodes  due  to  equilibrium  considerations  [11].  The 
total  potential  energy  for  the  overall  structure  is; 

P  =  V2{u}^[K]{u}  -  +  S  (1.66) 

In  the  current  study,  the  thermal  and  residual  effects  are  negli¬ 
gible  and,  therefore, 

=  S  =  0  (1.67) 

and 

P  =  l*{u}’'[Kllu}  -  {ul’^lLl^  (1.68) 

Comparing  this  relation  to  the  Griffith  criterion,  the 
total  potential  energy  P  corresponds  to  the  available  energy  for 
crack  growth.  The  first  term  on  the  right  side  is  the  elastic 
strain  energy  of  the  structure,  and  the  second  term  on  the  right 
side  is  the  work  done  by  external  forces.  Equation  (1.68)  can  be 
rewritten  as: 

-W  =  U  -  F  (1.69) 

where 

W  =  -P  (1 


U  =  l/2{u}^[K]{u} 


(1.71) 
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which  leads  to: 
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G 


[  {L}g  -  [K]{u}  ]  +  {u}’^ 


|T  3[K] 


(1.76) 


In  order  to  solve  equation  (1*80),  a  finite  element  model 
must  be  constructed  for  the  crack  problem.  Once  the  mesh  has  been 
created,  an  initial  solution  is  generated  with  a  crack  length,  a, 
and  an  applied  in-plane  normal  loading  corresponding  to  the  Mode  I 
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stress  state.  The  results  from  this  solution  are  the  nodal 
displacement  vector  {u}  and  the  overall  stiffness  matrix  tK]^  for 
the  applied  load  case.  A  second  solution  is  then  generated  in 
which  a  contour  of  nodal  points  encompassing  the  crack  tip  is 
selected,  and  they  along  with  the  elements  interior  to  the  contour 
are  rigidly  advanced  ahead  of  the  initial  crack  tip  by  a  small 
increment  Aa,  as  shown  in  Figure  4.  There  is  no  applied  loading 
in  this  solution,  and  the  result  is  the  overall  stiffness  matrix 

[^la+Aa  for  the  advanced  crack  case.  With  the  compiled  output 
from  the  finite  element  solutions,  the  stress  intensity  factor  can 
be  determined  by  first  making  the  approximation: 


3[K]  ^  tK]a+Aa  ~  ^  A[K] 

3  a  Aa  Aa 


(1-81) 


Then,  from  equation  (1.80): 

-E 


K-  =  [ 

2(l-v2) 


(1.82) 


It  is  important  to  note  that  due  to  the  symmetry  of  the 
crack  configuration,  only  the  upper  half  of  the  crack  is  actually 
modeled,  as  shown  in  Figure  4.  Therefore,  the  value  of  the  strain 
energy  release  rate  G  obtained  from  such  a  finite  element  solution 
is  half  of  that  expressed  in  equation  (1.79).  The  latter  value  is 
used  to  calculate  Kj  in  equation  (1.82)  for  the  overall  crack  con¬ 
figuration. 
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rj)  Contour:  ABCD 
rj  Contour:  EFGH 
Crock  Tip:  j 


X 


Figure  4.  Finite  Element  Model  Crack  Tip  Advancement 
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Equation  (1.82)  essentially  provides  the  same  result  as 
that  discussed  by  Parks  except  for  a  minor  variation.  In  the 
solution  by  Parks,  he  makes  the  following  approximation  [7]: 

3[K]  =  y,  aiki  ]  =  y  ^  (1.83) 

3a  i=l  8a  Aa 

The  reason  for  this  version  is  that  when  the  contour  Fq  is 
advanced,  the  only  change  in  stiffness  is  that  of  the  elements 
lying  between  the  contour  Fq  and  the  next  larger  nodal  contour  Fi 
which  is  also  shown  in  Figure  4.  Therefore,  equation  (1.83)  is 
the  summation  of  the  stiffness  matrix  variations  of  the  elements 
within  contours  Fq  and  Fi  which  are  the  only  elements  that  are 
deformed  in  the  process. 

The  expression  in  equation  (1.81)  was  used  in  this  study 
instead  of  that  of  Parks  in  equation  (1.83)  so  as  to  increase  com¬ 
putational  efficiency.  The  drawback  of  the  latter  is  that  in 
order  to  solve  equation  (1.83)  the  specified  element  stiffness 
matrices  must  be  separated  out  from  the  overall  matrices,  the  dif¬ 
ferences  computed,  and  then  these  element  matrix  variations  incor¬ 
porated  back  into  the  overall  matrix  so  that  the  rest  of  the 
matrix  operations  can  be  carried  out.  By  using  equation  (1.81) 
instead,  8[K]/3a  is  calculated  directly  from  the  finite  element 
output  for  the  overall  stiffness  matrices,  thereby  simplifying  the 
computations  considerably. 


CHAPTER  2;  FINITE  ELEMENT  MODEL  DEVELOPMENT 


2.1  Computer  Program 

The  development  and  implementation  of  the  finite  element 
models  for  the  crack  configurations  examined  in  this  study, 
require  the  use  of  three  separate  program  packages.  The  overall 
process  generally  involves  mesh  generation,  finite  element  com¬ 
putation,  and  subsequent  calculation  of  the  stress  intensity  fac¬ 
tor.  All  computation  is  performed  on  a  Burroughs  B7700  computer. 

The  first  computational  step  is  that  of  finite  element 
mesh  generation  which  is  accomplished  utilizing  the  CADOT  computer 
program  developed  by  Quigley  [12].  The  program  is  a  computer- 
aided  design  package  capable  of  complex  finite  element  mesh  design 
and  generation  in  two  and  three  dimensions.  Such  a  programming 
capability  is  considered  essential  for  optimum  efficiency. 

Once  the  mesh  has  been  created  and  compiled,  it  is  then 
input  into  the  finite  element  program.  The  program  used  in  this 
research  effort  is  SAPV  which  was  developed  by  Bathe,  et  al.  [13]. 
The  elements  which  are  utilized  are  the  two-dimensional,  four- 
node,  plane  strain  quadrilateral  elements.  The  output  which  is  of 
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interest  here  consists  of  the  nodal  displacement  vectors  {u}  and 
the  overall  stiffness  matrices  [K] . 

Due  to  the  large  amount  of  data  output  in  the  stiffness 
matrices,  the  solution  of  equation  (1«82)  can  be  a  potentially 
unwieldy  one.  To  handle  this  problem,  the  program  MATRIXDP  (see 
Appendix  A)  was  created  by  the  author  to  carry  out  the  com¬ 
putations  using  minimum  computer  storage  space  and  thus  maximizing 
the  data  handling  capabilities.  In  order  to  do  this,  only  the 
upper  half  of  the  symmetric  banded  region  of  the  overall  stiffness 
matrix  is  read  in,  with  each  number  correlated  to  its  symmetric 
location  when  called  for  in  the  matrix  calculations.  The  result 
is  that  only  a  portion  of  the  large  overall  matrices  require  space 
in  an  array  rather  than  storing  the  entire  overall  matrices  in  a 
redundant  manner.  The  final  output  from  this  program  is  the  com¬ 
puted  stress  intensity  factor  from  equation  (1.82). 

It  is  important  to  note  that  all  phases  of  computation 
previously  outlined  are  performed  in  double  precision.  This  will 
be  discussed  further  in  the  results. 

2 .2  Mesh  Design  Criteria 

In  this  research  effort,  three  two-dimensional  crack  con¬ 
figurations  are  modeled  with  finite  elements.  They  are  the  double 
edge  notch  (Figure  5)  the  single  edge  notch  (Figure  6),  and  the 
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center  crack  (Figure  7).  In  all  cases  the  cracks  are  completely 
through  the  plate.  The  plates  are  of  unit  thickness  and  subject 
to  plane  strain  conditions. 

The  most  important  consideration  in  the  development  of  the 
finite  element  meshes  is  the  optimization  of  the  mesh.  To  accom¬ 
plish  this,  a  number  of  criteria  are  specified.  First,  conven¬ 
tional  elements  are  used  exclusively  throughout  the  finite  element 
mesh.  The  reason  for  such  element  selection  is  that,  due  to  the 
universal  availability  of  the  element  types  used  in  this  study, 
the  method  is  easily  incorporated  into  most  finite  element  pro¬ 
grams  currently  in  use.  This  is  in  contrast  to  methods  which  uti¬ 
lize  specialized  singularity  elements  and  are,  therefore,  not 
easily  transferred  to  other  programs.  In  addition,  the  ease  of 
application  of  the  conventional  elements  as  compared  to  the  more 
specialized  elements  makes  the  conventional  elements  a  more 
natural  choice.  Specifically,  the  two-dimensional,  four-node 
quadrilateral  element  is  used  in  this  work.  This  is  in  contrast 
to  the  model  developed  by  Parks  [7]  which  makes  use  of  triangular 
and  quadrilateral  elements.  The  exclusive  use  of  quadrilaterals 
is  chosen  here  for  reasons  of  simplicity  as  well  as  accuracy.  It 
is  easier  to  generate  a  mesh  of  all  quadrilaterials  than  one  with 
a  mix  of  triangles  and  quadrilaterals.  As  the  problems  become 
more  complex,  this  becomes  an  increasingly  important  considera¬ 
tion.  In  addition,  the  quadrilateral  elements  provide  a  somewhat 
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higher  degree  of  accuracy  than  the  triangular  elements,  and 
convergence  can  be  achieved  with  a  coarser  mesh  made  up  of 
quadrilaterals  [11,14,15] . 

Another  criterion  for  the  mesh  optimization  is  that  the 
mesh  should  be  made  as  coarse  as  possible  while  still  achieving 
the  desired  accuracy.  This  is  extremely  Important  in  view  of  the 
fact  that  while  a  solution  may  be  theoretically  correct,  if  it 
requires  excessive  development  and  computational  time,  then  it  is 
of  little  practical  use.  The  accuracy  which  is  sought  in  the  pre¬ 
sent  study  is  within  5  percent  of  the  analytical  solutions  for  a 
relatively  coarse  mesh. 

Additionally,  the  configuration  of  elements  should  be  such 
that  the  mesh  can  be  generated  by  computer  rather  than  manually 
specified.  Here  again,  a  program  such  as  CADOT  [12]  is  indispen- 
sible  in  the  complex  problems. 

The  accuracy  of  a  coarse  mesh  such  as  one  with  the  pre¬ 
viously  discussed  characteristics  has  been  found  to  be  very  sen¬ 
sitive  to  the  mesh  geometry.  Based  on  work  done  by  Guydish  and 
Fleming  [8],  it  can  be  shown  that  very  good  accuracy  can  be 
achieved  with  a  properly  constructed  coarse  mesh.  They  also  show 
that  the  accuracy  of  a  well-conceived  coarse  mesh  geometry  is 
superior  to  many  ill-conceived  fine  mesh  geometries.  Based  on 
these  results,  the  finite  element  mesh  geometries  of  this  current 
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research  were  constructed  within  the  following  guidelines  [8] : 

1)  The  smallest  element  is  at  the  crack  tip 
with  an  optimum  size  of  0.005a. 

2)  Element  size  should  be  increased  away  from  the 
crack  tip  by  means  of  a  geometric  progression 
for  at  least  one  crack  length  in  all  directions. 

3)  The  element  aspect  ratio  (Length/Width)  should 
not  be  more  than  50  for  individual  elements. 


2 .3  Mesh  Configurations 

The  first  crack  problem  to  be  modeled  is  the  double  edge 
notch.  Due  to  the  symmetry  of  the  problem,  it  is  only  necessary 
to  create  a  mesh  for  one-fourth  of  the  overall  plate.  This  region 
is  indicated  by  the  shaded  area  in  Figure  5.  To  account  for  this 
symmetry,  the  appropriate  boundary  conditions  are  imposed  along 
the  edges  as  indicated  in  the  figure.  The  study  of  this  par¬ 
ticular  crack  configuration  involved  the  creation  of  four  dif¬ 
ferent  meshes.  Their  designations  are  M-2-2,  M-2-3,  M-2-4,  and 
M-2-7.  Each  of  the  four  meshes  has  a  different  geometry  so  that 
the  effects  of  these  variations  on  the  solution  accuracy  can  be 


ascertained . 
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Mesh  M-2-2  can  be  seen  in  Figure  8.  It  has  384  degress  of 
freedom  and  is  composed  of  204  nodes  and  176  elements*  It  is  a 
relatively  coarse  mesh,  and  one  which  is  easily  generated.  The 
smallest  element  is  0.0254  cm  (0.01  inch)  on  each  side  and  is 
found  at  the  crack  tip.  This  size  is  equal  to  0.005a  as  specified 
in  section  2.2.  The  multiplication  factor  for  the  geometric 
progression  of  element  size  away  from  the  crack  tip  is  approxi¬ 
mately  3.5  on  side  3.  The  multiplier  for  sides  2  and  4  is 
approximately  1.5,  and  for  side  1  it  is  equal  to  1.0.  The  largest 
element  aspect  ratio  is  approximately  141.  This  exceeds  the  ideal 
maximum  of  50  which  was  previously  specified.  However,  the 
adverse  effect  on  accuracy  is  negligible  as  will  be  shown  in  the 
results.  Only  two  elements  have  ratios  this  high,  and  they  are 
away  from  the  crack  tip.  The  other  elements  do  fall  within  the 
specified  range. 

Mesh  M-2-3  is  shown  in  Figure  9.  It  has  954  degrees  of 
freedom  and  is  composed  of  493  nodes  and  448  elements .  The  mesh 
geometry  in  the  region  of  the  crack  is  shown  in  more  detail  in 
Figures  10  and  11.  It  is  a  more  refined  mesh  than  M-2-2,  but  one 
which  is  just  as  easily  generated.  The  main  distinction  between 
the  two  is  that  the  geometric  multiplier  for  side  3  has  been 
reduced  to  approximately  1*4  in  M-2-3.  This  results  in  the  finer 
mesh  as  illustrated,  as  well  as  reducing  the  largest  element 
aspect  ratio  to  a  value  of  approximately  57 . 
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Figure  8.  M-2-2  Mesh  Geometry 
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Mesh  M-2-4  is  a  more  distinct  contrast  to  the  previous  two 
configurations*  It  is  found  in  Figure  12  and  it  has  2085  degrees 
of  freedom,  along  with  1066  nodes  and  1000  elements.  The  largest 
difference  in  this  case  is  that  it  does  not  employ  a  geometric 
progression  for  increasing  the  element  size  away  from  the  notch. 
Instead,  the  element  size  remains  constant  over  extended  regions. 
For  the  area  around  the  crack  tip  extending  one  crack  length  in 
the  X  and  y  directions  the  elements  are  all  square  with  a  side 
dimension  of  0.254  cm  (0.1  inch).  In  the  area  outside  this  region 
the  elements  are  all  uniform  rectangles  with  width  in  the  x- 
direction  of  0.254  cm  (0.1  inch)  and  length  in  the  y-direction  of 
5.08  cm  (2.0  inch).  While  this  mesh  is  more  refined  overall  than 
the  previous  ones,  it  is  coarser  around  the  crack  tip.  The  ele¬ 
ments  adjacent  to  the  crack  tip  measure  0.254  cm  (0.1  inch)  on 
each  side,  equaling  0.05a  rather  than  the  optimum  0.005a.  The  net 
effect  of  this  versus  the  increased  overall  refinement  on  the 
solution  accuracy  is  discussed  in  the  results.  In  this  mesh,  the 
largest  element  aspect  ratio  is  only  equal  to  20. 

Mesh  M-2-7  was  constructed  in  an  effort  to  approximate  how 
coarse  the  mesh  could  be  made  without  exceeding  the  specified  5 
percent  error  limit.  This  mesh  is  found  in  Figure  13,  and  has  162 
degrees  of  freedom,  with  88  nodes  and  70  elements.  For  this  case, 
the  geometric  multiplier  is  approximately  3.5  for  sides  2,  3,  and 
4.  Side  1  again  has  a  multiplier  of  1.0.  The  largest  element 
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aspect  ratio  is  approximately  140  which  is  comparable  to  that  of 
mesh  M-2-2. 

The  single  edge  notch  problem  and  the  center  edge  crack 
problem  are  also  modeled  in  a  manner  similar  to  that  just 
discussed  for  the  double  edge  notch  problem.  The  single  edge 
notch  problem  is  symmetric  such  that  one-half  of  the  plate  is 
meshed  as  shown  in  Figure  6.  The  symmetry  of  the  center  crack 
problem  requires  that  only  one-fourth  of  the  plate  be  meshed  as 
indicated  in  Figure  7 •  The  same  previously  developed  meshes  are 
also  compatible  with  these  two  crack  geometries.  The  significant 
difference  is  that  for  each  of  the  three  crack  geometries,  the 
boundary  conditions  are  different.  The  boundary  conditions  for 
these  crack  geometries  are  indicated  in  Figures  5  through  7 • 


CHAPTER  3:  RESULTS 


3.1  Double  Edge  Notch 

3.1.1  Theoretical  Solution 

The  analytical  expression  for  the  Mode  I  stress  intensity 
factor  of  the  double  edge  notch  problem  is  expressed  as  [9]: 

Kj  =  a('iTa)^^  f(a/b)  (3.1) 

where  the  factor,  f(a/b),  is  a  function  of  the  crack  half-length, 
a,  and  the  plate  half-width,  b  (see  Figure  5).  A  number  of 
expressions  for  f(a/b)  have  been  developed  in  the  literature.  The 
one  chosen  for  this  study  can  be  found  in  Reference  [9]  and  is  the 
following: 

f(a/b)  =  1.12  +  0.203(a/b)  -  1.197(a/b)^  +  1.930(a/b)^ 

(3.2) 

This  equation  was  derived  by  Brown  based  on  a  least  square  fitting 

to  Bowie's  results.  The  assumed  accuracy  is  within  2.0  percent 

for  a/b  ratio  <  0.7.  It  should  be  noted  at  this  point  that  the 

theoretical  and  numerical  results  for  Kj  in  this  chapter  are  non- 

Vl 

dimensionalized  in  the  form  K^/aa  purposes  of  comparison 

among  varying  load  cases  or  crack  lengths. 
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3.1.2  Effect  of  Crack  Advance  Increment  on  Solution  Accuracy 


A  significant  variable  in  the  stiffness  derivative  model 
is  the  increment  by  which  the  crack  is  advanced,  Aa.  Comparison 
of  the  approximation  made  in  equation  (1.81): 


3[K]  „  tK]a+Aa  ~  ^^^a  ^  A[K] 

3  a  Aa  Aa 

with  the  definition  of  the  derivative: 


(1.81) 


9f(c)  _  lim  f(c+h)  -  f(c) 

3x  h->0  h  ^  ^ 

indicates  that  the  approximation  of  equation  (1.81)  should  con¬ 
verge  to  the  exact  value  as  Aa  approaches  zero  in  the  limit.  The 

actual  application  in  this  problem  involves  specifying  a  suf¬ 
ficiently  small  value  for  Aa  in  the  finite  element  mesh  to  provide 
reasonable  accuracy.  This  was  investigated  by  generating  a  number 
of  computer  solutions  of  the  double  edge  notch  problem  with  dif¬ 
ferent  values  of  Aa  for  each  one.  For  all  cases,  the  mesh  is 
M-2-2,  the  loading  is  the  same,  and  the  inner  contour  size  r/a  is 
kept  constant.  In  the  expression  r/a,  the  term  r  is  the  approxi¬ 
mate  radius  around  the  crack  tip  of  the  inner  contour.  The 
results  of  these  computations  are  tabulated  in  Table  1,  and  pre¬ 
sented  graphically  in  Figure  14.  With  Aa  equal  to  0.254  cm  (0.1 


inch),  which  is  the  largest  increment  of  crack  advance  which  was 
tested,  the  solution  is  imaginary.  This  increment  produces  a 
change  in  the  stiffness  matrix  such  that  the  expression  within  the 
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Table  1.  Variation  of  Crack  Advance  for  Double  Edge  Notch 
Mesh  M-2-2,  r/a  =  0.254  cm  (0.1  in.),  a  =  75.8  MPa  (11  Ksi) 


[  ] 

^  172  "^theoretical 

oa 


2.06225 


Aa  (in.) 

K^/aa'/2 

Percent 

Error 

IxlO"^ 

2.72515x10"^  i 

— 

1x10"^ 

1.91393 

7.192 

1x10"^ 

2.01165 

2.454 

1x10“*^ 

2.02119 

1.991 

1x10"^ 

2.02214 

1.945 

1x10"^ 

2.02224 

1.940 

1x10“^ 

2.02225 

1.940 
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square  root  In  equation  (1.82)  is  negative  and  the  solution  con¬ 
sequently  imaginary  and  of  no  use.  Decreasing  this  increment  by  a 
factor  of  ten  results  in  an  error  of  7.192  percent.  By  continuing 
to  decrease  the  increment  value,  the  error  continues  to  decrease 
correspondingly  and  approach  a  value  around  1.940  percent  for  an 
increment  of  2.54x10“^  cm  (1x10  ^  inch).  At  this  point,  smaller 
increments  of  crack  tip  advancement  appear  to  have  insignificant 
effect  upon  the  accuracy.  Therefore,  all  further  computer  solu- 
tions  were  produced  with  2.54x10  cm  (1x10  inch)  selected  as 
the  optimum  value  for  Aa. 

3.1.3  Effect  of  Crack  Tip  Element  Contour  Size  on  Solution 

Accuracy 

A  second  factor  which  has  an  effect  on  the  solution 
accuracy  is  the  size  of  the  inner  contour  of  elements  To  which  is 
specified  around  the  crack  tip.  This  effect  was  examined  by 
generating  another  series  of  computer  solutions  of  the  double  edge 
notch  problem,  again  using  mesh  M-2-2.  For  this  series  of  com¬ 
putations,  the  crack  advance  increment  and  the  loading  are  kept 
constant.  The  results  are  shown  in  Table  2.  The  contour  size 
ratios  range  from  0.  (i.e.  the  single  node  at  the  crack  tip)  to 
0.5  with  a  corresponding  decrease  in  error  from  2.891  percent  to 
1.553  percent  as  illustrated  in  Figure  15. 


For  the  mesh  configuration  of  M-2-2,  0.5  is  the  maximum 
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Table  2 
Mesh  M-2-2 


Mesh  M-2-3 


.  Variation  of 

To  Contour  Size 

for  Double  Edge  Notch 

,  Aa  =  2.54x10“^ 

cm  (1x10  ^  in. 

),  a  =  75.8  MPa  (11  Ksi) 

r  I  1 

=  2 

.06225 

^  theoretical 

aa 

via. 

K  /aa^ 

Percent 

I 

Error 

0. 

2.00264 

2.891 

0.1 

2.02225 

1.940 

0.3 

2.02621 

1.748 

0.5 

2.03022 

1.553 

,  Aa 


2.54x10"'^  cm  (lxl0“^ 

[  ] 

'■  ]72'^  theoretical 

aa 


in. ) >  ®  “ 
=  2.06225 


72.4  MPa  (10.5  Ksi) 


0.5  2.06131 


0.046 


0.75 


2.06223 


0.001 


Percent  Error 
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possible  size  ratio  for  the  inner  element  contour.  The  more 

refined  mesh  M-2-3  can  have  a  contour  size  ratio  up  to  0.75.  The 

results  for  mesh  M-2-3  with  inner  element  contour  size  ratios  of 

0.5  and  0.75  are  presented  in  Table  2,  with  error  values  of  0.046 
percent  and  0.001  percent  respectively.  Again,  the  same  trend  is 

indicated  as  it  is  with  mesh  M-2-2.  The  results  indicate  that  as 

the  contour  size  To  increases,  the  solution  accuracy  increases 
correspondingly.  Intuitively,  this  is  a  logical  trend.  The 
region  close  to  the  crack  tip  will  show  the  most  significant  error 
between  the  theoretical  and  numerical  predictions  of  the  stress 
state  due  to  the  presence  of  the  stress  singularity  at  the  crack 
tip.  It  is,  therefore,  reasonable  to  assume  that  the  inner  ele¬ 
ment  contour  should  be  of  adequate  size  so  as  to  encompass  this 
region  within  the  rigidly  translated  elements,  with  the  deformed 
elements  then  at  a  sufficient  distance  such  that  the  singularity 
effect  is  minimized.  This  assumption  is  supported  by  the  results. 

3.1.4  Effect  of  Mesh  Refinement  on  Solution  Accuracy 

A  major  factor  which  influences  the  accuracy  of  the 
numerical  solution  is  the  mesh  geometry.  As  was  discussed  in 
section  2.2,  much  better  accuracy  can  be  achieved  with  a  well- 
designed  coarse  mesh  than  with  a  fine  mesh  that  has  an  inefficient 
geometry.  In  an  effort  to  demonstrate  such  effects,  a  series  of 
computer  solutions  were  generated  for  the  double  edge  notch 
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problem  for  a  variety  of  mesh  configurations*  The  results  are 
tabulated  in  Table  3*  As  is  indicated,  the  crack  advance  incre¬ 
ment,  and  the  element  contour  size  are  constant  throughout.  The 
applied  stresses  are  varied  somewhat  from  case  to  case  but  this 
has  no  effect  on  the  nondimensionalized  solution  or  its  accuracy. 

Comparison  of  meshes  M-2-2  and  M-2-3,  which  were  discussed 
in  section  2.3,  reveals  that  the  latter  mesh  has  essentially  the 
same  geometric  configuration  as  the  former.  The  distinction  bet¬ 
ween  them  is  that  M-2-3  has  approximately  2.5  times  the  number  of 
elements  as  does  M-2-2.  The  result  is  a  decrease  in  error  from 
1.553  percent  to  0.046  percent. 

Utilizing  the  same  basic  geometric  characteristics  of  mesh 
m-2-2  in  a  coarser  configuration  resulted  in  mesh  M-2-7 .  The 
results  are  that  mesh  M-2-2  has  approximately  2.5  times  the  number 
of  elements  that  M-2-7  has,  with  an  error  of  1.553  percent  for 
M-2-2,  and  an  error  of  5.549  percent  for  the  coarse  mesh  M-2-7. 

The  results  for  these  three  meshes  with  the  same  basic  geometric 
configuration  are  illustrated  in  Figure  16. 

The  results  indicate  that  an  increase  in  accuracy  can  be 
achieved  with  mesh  refinement  if  it  is  handled  properly.  The 
judgement  to  be  made  in  application  then  becomes  a  question  of  how 
accurate  a  solution  is  required  weighed  against  the  additional 
cost  of  using  the  more  refined  mesh.  In  this  specific  situation, 
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Table  3.  Mesh  Refinement  for  Double  Edge  Notch 


Aa  = 

=  2.54x10' 

^  cm  (1x10  ^ 

in.),  r/a  =  0.5 

^theoretical 

=  2.06225 

K  /oa^2 

Degrees  of 

Percent 

Mesh 

a 

Freedom 

Error 

(MPa) 

(Ksi) 

i. 

M-2-2 

75.8 

11 

2.03022 

384 

1.553 

M-2-3 

72. A 

10.5 

2.06131 

954 

0.046 

M-2-4 

69 

10 

1.86797 

2085 

9.421 

M-2-7 

69 

10 

1.94782 

162 

5.549 

Percent  Error 


5 


Figure  16.  Percent  Solutio 
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mesh  M-2-2  Is  sufficiently  accurate  in  view  of  the  minor  addi¬ 
tional  gain  that  is  made  at  the  expense  of  more  than  doubling  the 
mesh  size.  By  contrast,  comparison  of  the  results  for  mesh  M-2-4 
with  the  other  three  cases  exemplifies  the  problems  encountered 
with  poorly  designed  configurations.  Despite  the  fact  that  M-2-4 
has  1000  elements,  which  is  approximately  5.6  times  the  number  of 
M-2-2,  the  error  went  from  1.553  percent  for  M-2-2  to  9.421  per¬ 
cent  for  M-2-4.  If  M-2-4  is  examined  against  the  guidelines  for 
mesh  optimization  in  section  2.2,  it  can  be  seen  that  the  guideli¬ 
nes  are  not  adhered  to.  Specifically  the  size  of  the  elements 
around  the  crack  tip  are  too  large,  and  instead  of  a  geometric 
progression  in  element  size  growth  away  from  the  crack  tip,  the 
elements  remain  uniform  in  size.  The  result  is  that  despite  the 
much  denser  mesh,  the  accuracy  becomes  much  worse  due  to  the  poor 
layout,  and  computational  resources  are  inefficiently  utilized. 


3 .2  Single  Edge  Notch 

3.2.1  Theoretical  Solution 

For  the  single  edge  notch  problem,  the  stress  intensity 
factor  is  expressed  as  in  equation  (3.1).  The  f(a/b)  factor  which 
was  chosen  is  [9]: 
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f(a/b) 


0.265(1  -  a/b)‘*  + 


0.857  4-  0.265(a/b) 
(1  -  a/b)^''^ 


(3.4) 


which  was  derived  by  Tada.  The  assumed  accuracy  is  0*5  percent 
for  a/b  >  0.2. 


3.2.2  Numerical  Solution 


The  numerical  solution  of  the  single  edge  notch  problem 
was  produced  utilizing  the  optimum  parameters  which  were  deter¬ 
mined  in  section  3.1.  Mesh  M-2-3  is  used  to  model  the  crack  with 

~7  —7 

a  crack  advance  increment  of  2.54x10  cm  (1x10  inch)  and  an 
inner  element  contour  size  ratio  of  0.5.  The  results  are  tabu¬ 
lated  in  Table  4  with  a  resulting  error  of  0.858  percent. 

3 .3  Center  Crack 

3.3.1  Theoretical  Solution 

The  stress  intensity  factor  for  the  center  crack  problem 
is  also  expressed  as  in  equation  (3.1).  The  f(a/b)  which  was  cho¬ 
sen  is  the  following  [9]: 

f(a/b)  =  [(2b/iTa)  tan(Tra/2b) ] ^  (3.5) 

This  equation  was  derived  by  Irwin  and  is  based  on  an  approxima¬ 
tion  using  a  periodic  crack  solution.  The  assumed  accuracy  is 
within  5  percent  for  a/b  <0.5. 
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Table  4.  Crack  Configuration  Results 


Mesh  M-2-3,  Aa 

=  a.SAxlO""^  , 

cm  (1x10  ^  in.),  r/a 

=  0.5, 

a  =  72.4  MPa 

(10.5  Ksi) 

Crack  Configuration 

V2 

Kj/oa  ^ 

/  oa 

Percent 

Theoretical 

Numerical 

Error 

Double  edge  notch 

2.06100 

2.06131 

0.046 

Single  edge  notch 

4.98997 

4.94716 

0.858 

Center  crack 

2.00000 

1.95959 

2.021 
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3.3.2  Numerical  Solution 

The  numerical  solution  for  the  center  crack  problem  was 
also  arrived  at  by  employing  the  optimum  parameters  from  section 
3.1.  The  selected  mesh  configuration  is,  therefore,  M-2-3  with  a 
crack  advance  increment  of  2.54x10  ^  cm  (1x10  inch)  and  an  inner 
element  contour  size  ratio  of  0.5.  The  error  for  this  solution  is 
2.021  percent.  The  results  are  shown  in  Table  4. 

3  ,4  Sample  Computation 

A  complete  output  data  set  from  the  numerical  solution  of 
one  of  the  double  edge  notch  computations  has  been  reproduced  in 
Appendix  B.  The  data  corresponds  to  the  solution  generated  using 
mesh  M-2-7.  The  input  parameters  as  well  as  the  results  are  tabu¬ 
lated  in  Table  3.  Specifically,  the  data  is  the  output  from  the 
finite  element  program,  and  consists  of  the  overall  stiffness 
matrix  for  the  applied  loading  configuration,  the  corresponding 
nodal  displacements,  and  the  overall  stiffness  matrix  for  the 
advanced  crack  configuration. 

The  stiffness  matrix  ouput  includes  only  the  upper  portion 
of  the  symmetric  banded  region  of  the  actual  matrix.  From  this 
data,  the  program  MATRIXDP  generates  the  complete  set  of  terms  for 
the  matrix  calculations  as  described  in  section  2.1.  In  this 
solution,  the  stiffness  matrix  half -band width  including  the 
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diagonal  term  is  20.  The  stiffness  data  in  Appendix  B  is  read  by 
the  program  from  left  to  right  moving  down  the  page.  Reading  the 
data  in  this  fashion,  every  20  terms  makes  up  one  row  of  the 
matrix  half-bandwidth.  The  first  of  each  20  terms  corresponds  to 
the  diagonal  position  for  the  specified  row,  and  the  subsequent 
terms  correspond  to  consecutive  positions  to  the  right  of  the 
diagonal  term  in  the  array. 

The  nodal  displacement  data  is  arranged  with  three  rows  of 
data  for  each  nonstationary  nodal  point.  The  first  row,  reading 
from  left  to  right,  is  comprised  of  the  node  number,  the  x- 
direction  translation,  and  the  y-direction  translation.  The 
second  row  contains  the  z-direction  translation,  and  the  x- 
direction  rotation.  The  third  row  is  made  up  of  the  y-direction 
rotation  and  the  z-direction  rotation.  Since  the  mesh  was 
generated  in  the  y-z  plane,  only  the  y  and  z  translations  have 
values,  and  are  the  only  components  of  this  data  used  by  MATRIXDP 
to  generate  the  nodal  displacement  vector. 

The  data  in  Appendix  B  is  tabulated  as  it  is  read  in  by 
the  program  MATRIXDP.  The  first  data  file  which  is  read  by  the 
program,  consists  of  the  stiffness  matrix  data  for  the  applied 
loading  case,  on  pages  71  through  95,  and  the  data  for  the 
resulting  nodal  displacements  as  found  on  pages  96  through  99. 

The  stiffness  matrix  data  for  the  advanced  crack  case  is  on  pages 
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100  through  124,  and  comprises  the  second  data  file  which  is  read 
by  the  program. 

The  MATRIXDP  program  is  a  user-interactive  program.  For 
the  data  files  in  Appendix  B,  the  user  entries  are: 

Dimension  of  Matrix:  162 

Half  Bandwidth:  20 

Number  of  Nodes  (nonstationary):  87 

The  data  is  read  by  the  program  from  stored  pack  files  in  the  com¬ 
puter  system. 

An  important  factor  in  the  numerical  solutions  of  the 
problems  in  this  study  is  that  it  is  essential  that  all  com¬ 
putations  are  performed  in  double  precision.  Double  precision  in 
this  study  is  23  significant  figures,  as  shown  in  Appendix  B. 

This  degree  of  accuracy  is  necessary  in  order  to  generate  meaning¬ 


ful  results. 


CHAPTER  4 :  CONCLUSIONS 


The  stiffness  derivative  finite  element  technique  [7]  has 
been  utilized  in  this  study  to  determine  the  Mode  I  stress  inten¬ 
sity  factors  for  three  crack  configurations.  The  geometries  exa¬ 
mined  included  the  double  edge  notch,  the  single  edge  notch,  and 
the  center  crack.  The  primary  guidelines  in  developing  the  finite 
element  models  were: 

1)  To  use  standard,  universally  available  elements 
throughout  the  finite  element  mesh  including  those 
at  the  crack  tip. 

2)  Optimize  the  finite  element  mesh  geometry  for 
automatic  generation  by  the  computer. 

3)  Construct  the  finite  element  mesh  as  coarse  as 
possible  while  maintaining  the  required  accuracy. 

4)  Determine  the  computational  requirements  of  the 
"Stiffness  Derivative  Method,”  and  verify  the 
method's  accuracy  compared  to  established  analyti¬ 
cal  solutions. 
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The  results  indicate  that,  following  the  specified  guide¬ 
lines  using  the  "Stiffness  Derivative  Method,"  it  is  possible  to 
achieve  a  high  degree  of  accuracy  with  an  optimized,  relatively 
coarse  finite  element  mesh  composed  of  standard,  four-node,  plane 
strain,  quadrilateral  elements. 

The  design  of  the  mesh  configuration  was  found  to  have 
considerable  influence  on  the  accuracy  of  the  solutions  generated. 
For  a  well  conceived  mesh  design  with  degrees  of  freedom  ranging 
from  162  to  954,  the  error  ranged  from  5.549  percent  to  0.001  per¬ 
cent  respectively  for  the  double  edge  notch.  In  contrast,  a  mesh 
designed  with  less  optimization  of  the  geometric  configuration  and 
2085  degrees  of  freedom,  had  an  error  of  9.421  percent.  The 
results  for  this  mesh  design  resulted  in  a  waste  of  computing 
resources  and  insufficient  solution  accuracy.  Therefore,  con¬ 
siderable  savings  in  computing  time  along  with  a  high  degree  of 
accuracy  can  be  realized  with  a  well  developed  mesh  design. 

In  addition,  certain  parameters  involved  in  the  applica¬ 
tion  of  the  "Stiffness  Derivative  Method”  were  also  optimized  to 
maximize  the  solution  accuracy.  One  important  parameter  is  the 
size  of  the  increment  of  crack  advance  which  is  specified.  The 
trend  that  was  established  is  that  the  Increment  should  be  mini¬ 
mized  in  order  to  increase  the  solution  accuracy.  This  was 
demonstrated  for  the  double  edge  notch  problem  in  which  the 


63 


—2  —7 

increment  was  decreased  from  10  in.  to  10  in.  with  a  corres¬ 
ponding  decrease  in  error  of  7.192  percent  to  1.940  percent 
respectively. 

The  second  parameter  examined  for  its  effect  on  the  solu¬ 
tion  accuracy  was  the  size  of  the  element  contour  To  encompassing 
the  crack  tip.  The  results  indicated  that  increasing  the  radius 
of  this  contour,  increases  the  solution  accuracy.  For  the  double 
edge  notch  problem,  an  increase  in  the  element  contour  size  ratio, 
r/a,  from  0  to  0.5  reduced  the  error  from  2.891  percent  to  1.553 
percent. 


In  summary,  the  "Stiffness  Derivative  Method"  has  been 
employed  in  the  finite  element  models  of  selected  crack  con¬ 
figurations,  and  has  been  found  to  be  fully  compatible  with  the 
standard  four-node  quadrilateral  element  which  is  typical  of  most 
finite  element  programs.  The  mesh  design  which  produced  the 
highest  accuracy  was  found  to  be  relatively  easy  to  generate,  and 
utilizing  the  same  basic  geometry,  the  mesh  can  be  refined  to  pro¬ 
duce  the  desired  degree  of  solution  accuracy  within  reasonable 
computational  time  limitations.  The  results  of  this  study  are 
within  acceptable  error  limits  with  a  relatively  coarse  mesh.  In 
addition,  the  implementation  of  the  method  is  significantly 
simpler  than  other  methods  which  require  specialized  singularity 
elements,  extrapolation  of  data  from  numerous  generated  solutions, 
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or  changing  boundary  conditions  at  the  crack  tip. 

Based  on  these  favorable  conclusions  for  the  two- 
dimensional  crack  configuration,  the  next  logical  step  is  that  of 
a  three-dimensional  crack  configuration  such  as  the  penny-crack. 
The  application  should  differ  in  certain  respects  such  as  the  ele¬ 
ment  type  and  the  method  of  calculation  of  the  stress  intensity 
factor,  which  will  vary  along  the  crack  front.  The  three- 
dimensional  geometry  is  also  a  more  complex  problem  com¬ 
putationally.  For  this  reason,  the  results  of  this  study  related 
to  the  optimization  of  the  method  and  computational  efficiency, 
are  significant  in  the  development  of  an  accurate  and  com¬ 
putationally  feasible  three-dimensional  model. 
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$RESET  FREE 

$SET  SUPRS  LIMIT=100  LINEINFO 

FILE  2(KIND=PACK,MAXRECSIZE=14.BL0CKSIZE=420) 

FILE  3(KIND=PACK,MAXRECSIZE=23,TITLE='MATOUT' .PROTECTION=SAVE ) 

FILE  5(KIND=RSMOTE) 

FILE  6(KIND«REM0TE) 

DOUBLE  PRECISION  DELK(60000) , AK(60000) . BK(60000) . U(2500) 

DOUBLE  PRECISION  PROD1 (2500) . SIGL 

DOUBLE  PRECISION  XTRAN , YTRAN , ZTRAN , XROT , YROT , ZROT . DPOT , CPOT , TDPOT 
DIMENSION  ANS(12) 

INTEGER  A , B , C , D , S . HBAND ,E,F.G,H.P,0 
LOGICAL  L, ISCS 
DATA  ID0T/6H.  / 

DPOT=O.ODO 

E=0 

IBINFL=2 

ISCS=.TRUE. 

WRITE(6,2000) 

10  WRITE(6,2010) 

READ(5.2020)  ( ANS( I ) . I  =  1 . 11  ) 

ANS( 12)=ID0T 
CHANGE(IBINFL.TITLE=ANS) 

INOUIRE(IBINFL,PRESENT=L) 

IF  (L)  GO  TO  20 
WRITE(S.2030)  ( ANS (I ) , I = 1  ,  1 1 ) 

GO  TO  10 
20  E=E+1 
P=1 
0=2 
KK=0 

WRITE(S.2045) 

2045  FORMAT(/."  <ENTER>  DIMENSION  OF  MATRIX 11 X , 

+  "  HALF  BANDWIDTH",/, 1 IX, 

+  "  NUMBER  OF  NODES  (NONSTATIONARY)") 

READ(5./)  M, HBAND, NODEN 

N*HBAND*M 

F*HBAN0/2 

G=HBAND-(F»2) 

25  IF(F, LT. 1 )  GO  TO  33 
DO  30  H»1 , F 
IF(E. EO. 1 )  GO  TO  27 
READ(2,2050)  (BK( I ) , I =P , 0) 

READ(2,2055)  DBK,DDBK 
GO  TO  29 

27  READ(2,2050)  ( AK( I ) , I =P , 0 ) 

READ(2,2055)  DAK,DDAK 

29  P=P+2 
0=0+2 

30  CONTINUE 
33  0=P+G-1 

IF(P.GT.O)  GO  TO  37 
IF(E. EO. 1 )  GO  TO  35 
READ(2,2060)  ( BK ( I ) , I =P , 0 ) 

READ(2,20S5)  DBK,DOBK 
GO  TO  37 

35  READ(2,2060)  ( AK( I ) , I =P , 0 ) 

READ(2,2055)  DAK,DDAK 
37  P=0+1 
0=P+1 

IF  (I .LT.N)  GO  TO  25 
IF  (E.GE.2)  GO  TO  100 
DO  50  11=1 , NODEN 
READ(2,2070)  NNUM, XTRAN, YTRAN 
READ(2,2072)  ZTRAN, XROT 
READ(2,2074)  YROT, ZROT 
READ(2,2055)  DUM1 ,DUM2 
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IF  (YTRAN.EO.O. )  60  TO  40 

KK=KK+1 

U(KK)=YTRAN 

40  IF(2TRAN.E0.0. )  GO  TO  50 
KK»KK+1 
U(KK)«2TRAN 
50  CONTINUE 
CLOSE  2 
GO  TO  10 
100  CLOSE  2 

DO  500  S“1.N 
DELK(S)=BK(S)-AK(S) 

500  CONTINUE 

DO  700  I»1.M 
WRITE(3,510)  I 
510  FORMAT( "  I=" , 15) 

D=I+HBAND-1 

DO  600  U»I.O 

A=( (1-1 )«HBAND+d-I+1  ) 

IF  (d.GT.I)  GO  TO  580 
PR0D1(I)=DELK(A)*U(d) 

525  B=I-1 

IF  (B.LT. 1)  GO  TO  600 

C=I-HBAND+1 

IF  (C.LT.  1)  C  =  1 

DO  570  ISYM«C,B 

ASYM=( (ISYM-1 )»HBAND+I-ISYM+1 ) 

PRODK I )=PR0D1(I)+DELK(ASYM)*U(ISYM) 

570  CONTINUE 
GO  TO  600 

580  PRODK I )=PR0D1(I )+DELK(A)-U(d) 

600  CONTINUE 

DP0T=DP0T+(U(I)*PR0D1(I ) ) 

CPOT=U(I)*PROD1(I) 

WRITE (3. 650)  CPOT 
650  FORMATC  CPOT=  "  .  1PE  1 6 . 7  ) 

700  CONTINUE 

TDP0T=(-OP0T/2.OO0) 

SIGL=>DS0RT(-DP0T»1 .507/(9. 10-1*1 .00-7) )/2.0D4 
WRITE (6, 2080)  TOPOT 
WRITE(6.2085)  SIGL 

2000  FORMAT(/,"  STIFFNESS  DATA  INPUT  ROUTINE") 

2010  FORMAT("  <ENTER>INPUT  FILENAME") 

2020  F0RMAT(11A6) 

2030  FORMAT(/"***ERROR:  FILE  IS  NOT  PRESENT  ON  PACK" 
•/"FILE  “.11A6) 

2050  FORMAT (2030. 23) 

2055  F0RMAT(2F10.O) 

2060  FORMAT (2030. 23) 

2070  FORMAT (F6. 0,2030. 23) 

2072  FORMAT (036. 23, 030. 23) 

2074  FORMAT (036. 23, 030. 23) 

2080  FORMAT (//,"  POTENTIAL  ENERGY  CHANGE  =  ",1PD30.23) 
2085  FORMAT(//,“  K/(SIG*SORTA)  =  ",1PD30.23) 

STOP 

END 
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.537611 5245095 1 36605853 1 D+09 
.537307683629459023621020+09 


.7221 402 17551501 20526907D+07 
.  146712355720420140812780+07 


.268894903144812379285270+09 
.269198744024867016249560+09 
.  188762333829435403088300+10 
.  188702152324767555101  170+10 


145739637380385854873370+07 
.721 167499211466919329670+07 
■  . 143672259664195243802510+07 


.944589330762453257215140+09 
.  945191 145809 1 3 1 737086530+09 
.6911 14251370772415159590+09 
.  153156850861683356193560+09 


.  144644978004229529741880+07 
.  721167499211466919329690+07 

.  556008143297853474643080+05 
.  147305189372602451546860+07 


.268719199499224220340280+09 
.  345207920099826560846300+09 
.7713843748023237585101 10+08 
.  242559320707866664497590+10 
.5383852726344942945731 10+09 


.71921 8598579992248237090+07 
.  8572 1 2082015975230873220+04 
.  722002962386579550470010+07 
.  143079426012012933068450+07 


.943508736429419701436420+09 
.  121366282054987937378390+10 
.  27034049531695647 1738620+09 
.  198275162597745981151930+09 
.  42829458050782 1 9002 1 3240+08 


'  .719218598579992248237070+07 

■  . 857212082015975230833680+04 
. 722002962386579550470010+07 

). 

■  . 100166113945384054104870+06 
.  146306964169909162357650+07 


.  766566260346072060435450+08 
.969379327192409173571450+08 
.225701603699141462506480+08 
.6913653663121311 1 1728980+09 
.  152323468193025702569310+09 


.718389818441914492537710+07 
.  176751467946954407469560+05 
.  723411376999087124584310+07 
.  144077651214706222257640+07 


.  268657382483558032263670+09 
.345327522633821969757270+09 
.  773267656348750520801800+08 
.  6101802233 1 5354991093600+08 
.  100673454849106106470740+08 


'  .718389818441914492537710+07 
.  176751467946954407468450+05 
.723411376999087124584310+07 
). 

■  .218729350460727771649720+06 
. 141862017585926764086530+07 


OOOOO  IIIOOOOO  IIIOOOOO  I  lOOOOOJ  IIIOOOOO  I  lOOOOO 


73 


0. 

c. 

0- 

0. 

0. 

“ . 9979986037708656 1476S43D+07 
- . 4620483265 1 53402752 1 536D+08 
18321008791764518998298D+08 
. 13521380417019281734946D+09 
0. 

0. 

0, 

0. 

0. 

0. 

- . 57314154540015679323496D+08 
- . 135031478639825194512330+09 
0. 

.397822072480211708571250+08 

0. 

0. 

0. 

0. 

0. 

0. 

- . 164212069483473646164590+08 
- . 382463797676357698 1 76560+08 
0. 

. 704655482929487995088150+09 
- . 703878786359643629604490+09 
0. 

0. 

0. 

0. 

0. 

0. 

. 825924796845358762833470+08 
- . 835448798999684007 1 1 9680+08 
. 2473229222183611 12138250+10 
- . 247247188641 183170280960+10 
0. 

0. 

0. 

0. 

0. 

0. 

. 291044777375216297526490+09 
- . 292882707480029271878120+09 
. 906873202668587034374570+09 
- .2002414869467754252741 10+09 
0. 

0. 

0. 

0. 

0. 

- .825706986 1 1 3896733059390+08 
. 105181447231598778183110+09 
- . 247162280227294229263450+08 
.317875977353872766644680+10 
- . 705488749302310215362390+09 
0. 

0. 

0. 

0. 

0. 

- .289462002123769309505780+09 


0. 

0. 

0. 

0. 

-.621 1396966503647 194004 10+07 
.499464464234809227182050+07 
. 89455005382596209 1 807930+07 
0. 

.315407659199589350628520+07 

0. 

0. 

0. 

0. 

0. 

0. 

- . 548041479914675079176670+07 
.547321815783965014456620+07 
0. 

0. 

0, 

0. 

0. 

0. 

0. 

- . 548041479914675079176670+07 
. 25693720039934962984 1 370+07 
0. 

0. 

- . 3936451549736 1 1578555490+05 
. 689092109250866679453320+04 
0. 

0. 

0. 

0. 

0. 

0. 

143364874658394104764330+07 
. 721185854674642072365630+07 
. 689092109250866679444990+04 
0. 

0. 

0. 

0. 

0. 

0. 

. 147019740726221279850950+07 
.721  185854674642072365640+07 
0. 

142978335471557338429320+06 
- . 100726 1 19658293832093560+05 
0. 

0. 

0. 

0. 

0. 

- . 719199383493481605185720+07 
. 832705323770635386872100+05 
. 723567364202726605240270+07 
- . 1007261 19658293832092660+05 
0. 

0. 

0. 

0. 

0. 

- .719199383493481605185720+07 
. 832705323770635386874810+05 
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. 373885998S2467722852264D+09 
-,850374265826832717254580+08 
.264376611 194903828667010+09 
- .54459198912551 1256381260+08 


.  2374 1 58797 1 03904682831 70+08 
.254897671073466694749420+08 
.914319753626519678573460+07 
.  907635540983472605665830+09 
.  199337953046133902205360+09 


.8161 659209384788938 1 1 590+08 
,  105386582573807342306950+09 
.  259244525 163302199978 1 20+08 
.  95709753303778 162877 1720+08 
.719803582211941031092270+07 


.  8 1 66954686 1 776 1 2206 1 5080+07 
.864689202186390137632800+07 
,91460191 4096808652006690+07 
.  263008530 1 7 1 059 1 05 1 48890+09 
.  537857372322032169527360+08 


.  2250316449753 1 5141726090+08 
.253535058607472252321650+08 
.  104135269899688487851480+08 
.966776991295953341977500+08 
.245740494820997861585980+08 


.816413155292974730073590+07 
.  549 1 4035053050538045 1 1 60+08 
.231039970846111806012990+08 
.  9497010758 1 1907184 127240+08 
.808001615601318866326240+07 


.699060185559092650317330+07 
.792719429257595517940030+07 
.  101906555061317031703670+08 
.  2834403725150132967 1 1260+09 
.  1 12151672790886777118280+09 


. 723567864202726605240270+07 
0. 

- .400820847410756141089510+06 
-. 102662531517888671845820+06 
0. 

0. 

0. 

0. 

0. 

-.716820148675552448921380+07 
.295491384093528735276360+06 
.731861 162873931854727490+07 
- . 102662531517888671845790+06 
0. 

0. 

0. 

0. 

0. 

- .716820148675552448921360+07 
. 29549 1 384093528735276400+06 
.731861162873931854727510+07 
0. 

126531546649374670987370+07 
- . 437370550355949278049940+06 
0. 

0. 

0. 

0. 

0. 

- . 7085365396 1 7471422059420+07 
. 104143332270833429098710+07 
. 760924324396408533505680+07 
- . 437370550355949278049900+06 
0. 

0. 

0. 

0. 

0. 

- . 7085365396 1 747 1 422059420+07 
. 104143332270833429098720+07 
, 760924324396408533505700+07 
0. 

- . 407519937107802552090220+07 
-. 151069264154157608954470+07 
0. 

0. 

0. 

0. 

0. 

- .679507428271995428576190+07 
. 3460406 1 3850205924403010+07 
.853122696531472885913220+07 
- . 151069264154157608954460+07 
0. 

0. 

0. 

0. 

0. 

- .679507428271995428576180+07 
. 3460406 1 3850205924402980+07 
.8531 2269653 1 4728859 1 32 1 0+07 
0. 

- . 1 4523025395089102 1 322490+08 
- . 5782664880975044 1 6406560+07 
0. 

0. 
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0, 

0. 

0. 

- . 22 1032775557 12400468090D+08 
-.202394457116635350991000-^09 
- ,808191671212963857982010+08 
. 101037249643131 163191240+09 
. 2894269203 1 6660874020030+08 
0. 

0. 

0. 

0. 

0. 

-  .  676223948268 1 7 1 886386820+07 
- . 565550426106160356699230+08 
- . 256342304867362068697340+08 

.497472446136152583856050+09 
. 140375599292062531563280+09 
0. 

0. 

0. 

0. 

0. 

- . 7975080 1 3 1 46920364208060+08 

-  .  326594443301215718636710+09 
- .840622567781758159460610+08 

. 150267156252728767414430+09 
. 38 1 40528588329205 1 762350+08 
0. 

0. 

0. 

0. 

0. 

- . 221861068151255414244120+08 
-.94131 6785999052572780610+08 

-  .  265592151808492755732930+08 
. 304147519032189171453370+09 

0. 

0. 

0. 

0. 

0. 

0. 

- . 830 1 295 12871 945656058 1 30+08 
- . 167943297188875178621440+09 
0, 

. 899807 1 6270664 1 619086250+08 
0. 

0. 

0. 

0. 

O. 

0. 

- .231 160713490872622330820+08 

-  .  484377849505058160358240+08 
0. 

.21452461521 1724451608750+09 
- .212459004604689727029470+09 
0. 

O. 

0. 

0. 

0. 

0. 

.2254898251504493351 17470+08 
- .246363741952258610684400+08 


0. 

0. 

0, 

- . 587422599873185901788960+07 
. 122470922436264830980390+08 
. 120369472887559776244400+08 
- . 578266488097504416406570+07 
0. 

0. 

0. 

0. 

0. 

- . 587422599873185901788950+07 
. 122470922436264830980400+08 
. 120369472887559776244390+08 
0. 

- .2478930358361251 14966440+08 
- .672692796794735018794710+07 
0. 

0. 

0. 

C. 

0. 

- .2372584583581 10153062250+07 
. 197120555330874555025020+08 
. 122335862633099479954870+08 

-  .  672692796794735018794750+07 
0. 

0. 

0. 

0. 

0. 

-  .  237258458358 1  10153062250+07 
. 197120555330874555025020+08 
.  122335862633099479954870+08 

0. 

- . 142019731379059143461 180+08 
0. 

0. 

0. 

0. 

0. 

0. 

-  .  21747864  1425830188428630+07 
.  1  15888149778584492018580+08 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

- .217478641425830188428640+07 
.868496882401229535570570+07 
0. 

0. 

- .69501 1675606573573931170+05 

-  .455024126448895550292660+05 
0. 

0. 

0. 

0. 

0. 

0. 

-. 138988105384574385522620+07 
. 722668 106 172389402099620+07 


OOOOt  IIIOOOOO  I  lOOOOOl  iiiOOOOOl  I  lOOOOOl  »  iOOOOOJ  I  oooooo 


76 


. 74882608468880650536984D+09 
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0. 
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0. 

0. 
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0. 
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0. 

0. 
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0. 

0. 
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0. 
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0. 
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0. 
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0. 
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0. 

0. 

0. 
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0. 
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0. 
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0. 

0. 

0. 

0. 
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,  2210327489 1 46 1 828980476D+08 
.202394454613073281147760+09 
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.  101037249160264775823700+09 
,  289426920316660874020030+08 


.676223838714316996393520+07 
.  565550421 274980067773700+08 
.  256342304867362068697340+08 
.497472446136152583856050+09 
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.  7975080 1 3 1 46920364208060+08 
.  326594443301215718636710+09 
.  84062256778 1 758 1 594606 1 0+08 
.  150267156252728767414430+09 
.  38140528588329205 1 762350+08 


.221861068151255414244120+08 
.  94 1 316785999052572780610+08 
.265592151808492755732930+08 
.  304147519032189171453370+09 


.  830 1 295 1 287 1 94  5656058 1 30+08 
.  167943297188875178621440+09 

.  899807162706641619086250+08 


.23116071 3490872622330820+08 
.  484377849505058 1 60358240+08 

.214524615211724451608750+09 
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. 225489825 1 50449335 1 1 7450+08 
- .246363741952258610684390+08 


0. 

0. 

0. 

- . 587422559240571654536680+07 
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0. 

0. 
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0. 

0. 

0. 

0. 
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0. 

0. 

0. 
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0. 
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0. 

0. 

0. 

0. 

0. 

0. 
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. 255436 1 24896865539205550+08 
. 124751087124327178573440+08 
- . 1034553854543253457028 10+08 
0. 

0. 

0. 

0. 

0. 

.864279030739663702250310+06 
. 255436 1 24896865539205560+08 
. 124751087124327178573440+08 
0. 

- . 198417311 101749724162160+08 
0. 
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.99969352774578926691427D+08 
.  18270573748325588332331D+09 


188884046733634419529250+07 
1  1  1575550356216078305970+08 


104200530040 1 5336334 1790+09 


.  280203830350265071023230+08 
.524740156481686264372640+08 

.656272089866760124276410+08 
.  583446440083283272440280+08 


.  188884046733634419529250+07 
.825370888177545398444410+07 


.288100162173630162948120+07 
.892952570883608272805610+07 
.214580279185877144932280+09 
.  212816640247662545386550+09 


.  2248620455 1931 10607 1 8050+08 
.  252000368729934019744280+08 
.  100260256631917973653730+09 
.  90328032265 1 857141945900+07 


- . 173110753985502972692790+06 
-.1091 33578759842547576540+06 
0. 

0. 

0. 

0. 

0. 

0. 

131807752092814998810850+07 
.726337294924352223713600+07 
- . 109133578759842547576270+06 
0. 

0. 

0. 

0. 

0. 

0. 

.  158576863291800385804450+07 
. 726337294924352223713580+07 
0. 

- .800705891 109322163173900+06 
- . 348017061024377390503870+06 


.838931857669547098182900+07 
.  834706847 1 654874801299 1 0+07 
.  878708677884475522499740+07 
.  28 1 303804409520 1 56432590+09 
.  590305040673801097509910+08 


- .714019219249916393789720+07 
.588800794726249031428180+06 
.741907594509337914412630+07 
- . 348017061024377390504070+06 
0. 

0. 

O. 


.  23276334 1 3 18250395396220+08 
.  256943303904437 162802980+08 
.974128503479094219892950+07 
.  97062 1 5 17316794460264490+08 
.231277103313184331740180+08 


- . 828325242561074145660010+07 


-.7140192 192499 1 6393789710+07 
. 588800794726249031428690+06 
. 741907594509337914412630+07 
0. 

- , 280552746548332062653500+07 
-.111 171959939847316907470+07 
0. 

0. 

0. 

0. 

0. 

- . 698325075906 1 20980827830+07 
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- . 548541088057C3001895560D+08 
- . 224465884074421026228980+08 
. 100827325098609082012960+09 
- . 982769 1 929760088 1 2097250+07 


. 782735852518729842382380+07 
.  748010642827339729194990+07 
.  924785042271716594824070+07 
.  240705 1 7258497 1 292490670+09 
.  829397529606067534 1 36300+08 


.  220594472414466957191000+08 
.  175345042843452361656380+09 
.  62803420237 1 52 1632716260+08 
.  9084430 1 346 1 243 1 8652 1 700+08 
.  187364963454405509057430+08 


.  736797082928560287758240+07 
.  477180392054678265234070+08 
.  195731365809382768172480+08 
.  589108244713652448704680+09 
.  191265349378995497255930+09 


.627257037543626699961390+08 
.  4 1 33255229 17104 153607600+09 
.  129946271622782872476860+09 
.  176674065514086954  1  19320+09 
.  530553609963806825553940+08 


.  177820347859820643520330+08 
.  1 1 7594489074209467 1 77590+09 
.  384969639020978574 1 55820+08 
.  1 19455809436713896281 170+10 
.  360868379682357546293790+09 


.  130391855261828209036680+09 
.  734513820786466270629320+09 
.  209107695602475280531070+09 
.  347094058977090209408110+09 
.  103008225866985128382570+09 


. 196332294748688461305010+07 
.792938061969261419802310+07 
-.111 171959939847316907470+07 
0. 

0. 

0. 

0. 

0. 

- . 698325075906 1 20980827830+07 
.  196332294748688461304990+07 
. 79293806 1 96926 1 4 1 98023 1 0+07 
0. 

- .900431280935357691718810+07 
-.33191 552982 1 939797002890+07 
0. 

0. 

0. 

0. 

0. 

-.646897926507694997050010+07 
. 600039355500026774020040+07 
. 93 14369 1 770385633 1 34 1 390+07 
- .331915529821939797002890+07 
0. 

0. 

0. 

0. 

0. 

- . 646897926507694997049990+07 
. 600039355500026774020040+07 
.931436917703856331341410+07 
0. 

- . 231957871 101690033893350+08 
-  .772816195763232018715440+07 
0. 

0. 

0. 

0. 

0. 

- . 50734704 1 700658469158920+07 
. 140741962333449370056670+08 
. 1 16002485946761299276820+08 
- .772816195763232018715430+07 
0. 

0. 

0. 

0. 

0. 

- . 50734704 1 700658469 1 58920+07 
. 140741962333449370056660+08 
. 116002485946761299276810+08 
0. 

-.492557289251020368855360+08 
- . 151907712334310389614270+08 
0. 

0. 

0. 

0. 

0. 

- . 276979951097714059804050+07 
.249882676037159901209240+08 
. 142927725345227940236520+08 
- . 151907712334310389614270+08 
0. 

0. 

0. 
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.  368560021077483 1 573 1 223D+08 
.  2102075 1 6802037 1  1 572589D+09 
.  610936065426490388044330+08 
.  104590332773934144733930+10 
.  205502998928787928975720+09 


.  210062885486868195402070+09 
.  563382530402223912116090+09 
.  106178550254006331602950+09 
.304141738363176763799100+09 
.57621 6669658658653 19191 0+08 


.  59598863958 1 1 78852 1 304 1 0+08 
.  161036453991795547298770+09 
.317492814674387902095440+08 
.  370744755539956997967800+09 


.  106323462022210930744520+09 
.  187144983693230792087730+09 

.  108886798391266333718170+09 


.  300209584377684439805330+08 
.  532036045718664737195770+08 

.  344492538 127707978 1 20440+08 
.884777494941628047647520+07 


.  1  14340391561878187353880-^08 
.727729268468858448174560+07 
.  655298278466434762049500+08 
.  586096678034068606033410+08 


. 72772926846885844817455D+C7 

0. 


0. 

0. 

- . 276979951097714059804060+07 
. 249882676037 1 59901209230+08 
.  142927725345227940236520+08 
0. 

- .446475224148703721503390+08 
- . 892840097124524708674760+07 
0. 

0. 

0. 

0. 

0. 

- . 599505484765124477159110+05 
. 191778043381618631184130+08 
. 107966273255998985693700+08 
- . 892840097 1 24524708674740+07 
0. 

0. 

0. 

0. 

0. 

- .599505484765124477159610+05 
. 191778043381618631 184140+08 
. 107966273255998985693700+08 
0. 

- . 160470887644577071 164860+08 
0. 

0. 

0. 

0. 

0. 

0. 

- . 358361721575295266615840+07 
. 783028935724773502760390+07 
0. 

0. 

0. 

0. 

0. 

0. 

0. 

- . 358361721575295266615840+07 
. 4926443203401581 18145100+07 
0. 

0. 

- . 275603497865342793429920+06 
-. 151246831629407482624040+06 
0. 

0. 

C. 

0. 

0. 

0. 

- .865912275220752924783340+07 
0. 

-  .  151246831629407482624020+06 
0. 

0. 

0. 

0. 

0. 

0. 

, 160782735823347972031410+07 
- . 905630018232821620406470+07 
0. 
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.971 13529383409365857180D+08 
.  23090466460555463243 1 590+08 


.848283101312374555417390+07 
.  669945210190998670736160+06 
.  745033462705919244023410+07 
.  100332698330587598222950+09 
.  99301 1151 124428630504320+07 


.  669945210190998670736260+06 
.745033462705919244023390+07 

.  239407835 148501561 946 1  1 0+09 
.  820010420051903823381030+08 


.  6952894948 1 34 1 5889798 1 00+07 
.203433092195053589911810+07 
.  791642558646265797341330+07 
.900828252566476620657100+08 
.  182162699025862321604290+08 


.  6952894948 1 34 1 5889798 1 1 0+07 
.  203433092 1 950535899 1 1 780+07 
.791642558646265797341350+07 

.601690773684571841983890+09 
.  195438529538936275438730+09 


.  64855 1031669935663904860+07 
.  490745940923274966926840+07 
.874377070011747952393760+07 
.  179040568550176442094500+09 
.538484199201501985566180+08 


64855103 1 669935663904820+07 
490745940923274966926850+07 
874377070011747952393770+07 

1 10965020146044262534840+10 
322535630759532324089390+09 


- . 1 2224502084 1 303506085 1 50+07 
“ .489073674419738001 169100+06 
0. 

0. 

0. 

0. 

0. 

- . 55008061782744 1054869240+08 
- . 223054814925789985968100+08 
0. 

- .489073674419738001 168710+06 
0. 

0. 

0. 

0. 

0. 

-.71 26432 1076965720266 1 1 60+07 
- . 747037958586193400443550+07 
- . 858556345738403907585070+07 
0. 

-.3951 57620137263264751090+07 
- . 14706313119899781 1254690+07 
0. 

0. 

0. 

0. 

- .221718885252650457580550+08 
-. 174451940838516492473950+09 
-  .621748714244714165057110+08 
0. 

- . 1470631311989978 11254700+07 
0. 

0. 

0. 

0. 

0. 

- . 794471618358744376298170+07 
- . 471315957755664244507650+08 
-. 187033093964857973388900+08 
0. 

- . 108607713702833740073550+08 
- . 369062086806381768089040+07 
0. 

0. 

0. 

0. 

-.621 5296 10820282399464650+08 
-.422965585114538009593670+09 
- .  1334944640268751  15995190+09 
0. 

- . 36906208680638 1 768089050+07 
0. 

0. 

0. 

0. 

0. 

- . 180943439452807356820290+08 
- . 1 19582123724647749721920+09 
- . 386808662888 1 73305839 1 50+08 
0. 

- . 220824499436813188437930+08 
- .667373817598964021692530+07 
0. 

0. 

0. 

0. 


Ill 


0. 

- . 56559747712886773381681D+07 
.8054356094881 10919370360+07 
. 942966 1 76473598516503290+07 
. 3202052 1 7668050942449 1 20+09 
.91071 380437276201 1 242510+08 
O. 

0. 

0. 

0. 

- .565597477128867733816790+07 
. 805435609488 1 10919370340+07 
.942966176473598516503290+07 
0. 

. 163369238914015232653430+10 
. 446088204696889219104300+09 
0. 

0. 

0. 

0. 

0. 

- .496814060412979827678580+07 
. 104097286550282176642780+08 
. 993903661495444769693550+07 
.46901820021 1 185363600820+09 
. 126747802278858732876170+09 
0. 

0. 

0. 

0. 

- . 4968 1 40604 1 2979827678560+07 
. 104097286550282176642780+08 
. 99390366 1 495444769693540+07 
0. 

. 1 15813814204381230298120+10 
. 2 1 3545 1 6430610806822 1 330+09 
O, 

0. 

0. 

0. 

0. 

- . 445774683788787250132090+07 
. 68887 133191 6282943996 1 00+07 
. 845 137326302656093460300+07 
. 3335998 1 5796703008069 1 50+09 
. 594869055779183378515310+08 
0. 

0. 

0. 

O. 

“ .4457746897887872501321 10+07 
.688871331916282943996080+07 
.845137326302656093460290+07 
0. 

. 376209 1 20559100035377950+09 
0. 

0. 

0. 

0. 

0. 

0. 

- . 594946375302045654710860+07 
.363955977031654922140290+07 
0. 

. 109572002440650323428560+09 
0. 


- ,  1336616336794539C788675D+09 
.  694756900506699093004820+09 
192685029156673698069020+09 

0. 

- . 6673738 17598964021692550+07 
0. 

0. 

0. 

0. 

0. 

- . 38 12580644768 1 9203945270+08 
197512155624242604730210+09 
- . 554634281906567572791890+08 
0. 

- . 3535182277 1 596504731 3680+08 
- .988362936818232973702310+07 
0. 

0. 

0. 

0. 

- . 193019929119247594189130+09 
-.898131 02075206969 1 267910+09 
- .236642296829139245762660+09 
0. 

-.988362936818232973702340+07 

0. 

0. 

0. 

0. 

0. 

- .549561800101077859586250+08 
- .255620186165527192634200+09 
- . 679576720894324606539050+08 
0. 

- .260122784635127679640780+08 
- .487339150953812464804830+07 
0. 

0. 

0. 

0. 

- .237065154459175735428300+09 
- . 593814190953362129249120+09 
-. 108078477607361602237740+09 
0. 

- . 487339 1 509538 1 2464804840+07 
0. 

0. 

0. 

0. 

0. 

- .674755019015522832791340+03 
- . 168700730230814909976880+09 
- .315072725488998554571060+08 
0. 

- . 853977503675944777706740+07 
0. 

0. 

0. 

0. 

0. 

- . 108150973332967940126550+09 
- . 188279777585003039782050+09 
0. 

0. 

0. 

0. 
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.  59494637530204565471088D+07 
.73571361647039537525019D+06 


.  343770701 13249548640835D+08 
.430211422042248159414140-15 


.848283101312374555419720+07 

.  12531340700914  1387709500+09 
.409080688859183933379010+08 


.  8659 1 2275220752924785640+07 
.669945210190998670752330+06 
.  6952894948 1 34 1 5889797300+07 
.685671686044426611 109370+08 
.  396232436944377391818440+07 


.  669945210190998670752380+06 
.  6952894948 1 34 1 5889797300+07 

.3538191 35599660300659580+09 
.  123093562340668768037800+09 


.  7450334627059 1 9244024 1 90+07 
.203433092195053589912620+07 
.  64855 1 03 1 669935663905240+07 
.  1  1  1275566872439791459410+09 
.  323851151  145788330399420+08 


.  745033462705919244024  190+07 
.  2034330921950535899 1 2650+07 
.  64855 1 03 1 669935663905280+07 

.  8479225 1 1 324970049477680+09 
.  26658810677 1 287850457730+09 


0. 

0. 

0. 

0. 

- . 309246097320693078 1 7 1 960+08 
-  .  531814122471940895337740+08 
0. 

0, 

.  54566702253735799384 1 380+07 
0. 

0. 

0. 

0. 

0. 

0. 

- . 114340391561878187354900+08 
. 7 1 26432 1076965720266 1 140+07 
0. 

. 162647673107585433172060-13 
. 7993605777301 1 2709105020-1 4 
0. 

0. 

0. 

O. 

0. 

- . 5500806 1 782744 1 054870390+08 
-.22171 888525265045758 1 1 40+08 
0. 

. 793809462606986926402900-14 
0. 

0. 

0. 

0. 

0. 

- . 727729268468858448 1 74590+07 
- .747037958586193400450650+07 
- . 79447 1 6 1 8358744376298420+07 
0. 

. 849320613838244753424080-14 
- .416333634234433702658850-14 
0. 

0. 

0. 

0. 

- . 2230548 1 4925789985968750+08 
-. 174451940838516492474570+09 
- . 62 1529610820282399466250+08 
0. 

-.41078251911 1307919956730-14 
0. 

0. 

0. 

0. 

0. 

- . 858556345738403907587240+07 
-.47131 5957755664244509460+08 
-, 180943439452807356820760+08 
0. 

. 32 1 964677 1 4 1 295396722850- 1 4 
. 466293670342565746977930-14 
O. 

0. 

0. 

0. 

-.62174871 42447 1 4 165058690+08 


- . 791642558646265797340940+07  - . 422965585 1 14538009594330+09 


OOOOil  O  OOOOOII  OOOOfi  O  OOOOOII  OOOOII  O  oooooiido 


113 


-.490745940923274966927200+07  133661633679453907886900+09 

.565597477128867733816330+07  0. 

24572591 1977213449311 120+09  .470457006684910084004500-14 

747114520458386525598820+08  0, 

0. 

0. 

0, 

0. 

791642558646265797340920+07  187033093964857973389310+08 

490745940923274966927180+07  1 195821237246477497221 10+09 

565597477128867733816340+07  - , 381258064476819203945700+08 

0. 

13906587981 1409414446040+10  - .478420596780777286349800+02 
385295683531169756223490+09  - . 196314246620655954433920+02 

0. 

0. 

0. 

O. 

- , 133494464026875115995340+09 
874377070011747952394240+07  - .694756900701098975946590+09 
805430827840415928473610+07  -. 193019930035755856010530+09 
496816026113718851114230+07  0. 

398666638291670499804940+09  -. 196314246620655954780870+02 
108920313651147562529800+09  0. 

0. 

0. 

0. 

0. 

874377070011747952394230+07  - . 386808662888173305839590+08 
805430827840415928473630+07  - . 197512153678239016097250-^09 
4968160261 137188511 14220+07  - . 549561820984326469406750+08 

0. 

17971 1068570552290021980+10  - . 306777387954257346400210+02 
473356038959273058579520+09  .832667268468867405317660-15 

0. 

0. 

0. 

0. 

- . 192685028430787131821740+09 
942964216085 1 55574609460+07  - . 898 1 31020938083723007470+09 
104096979497491895918860+08  - . 237065 1 54459 1 75735428390+03 
445774689788787250132030+07  0. 

514023662767676336275020+09  .915933995315754145849430-15 

134174251187896836693700+09  0. 

0. 

0. 

0. 

0. 

942964216085155574609460+07  - . 554634245 1 4492 1 452 1 86420+08 
104096979497491895918860+08  - . 2556201 84932940073808220+09 
445774689788787250132000+07  - . 674755019015522832791590+08 

0. 

118846231343807271131320+10  . 166533453693773481063550-14 

215396569833592860384140+09  .277555756156289135105910-15 

0. 

0. 

0. 

0. 

-  .236642296829139245762750+09 
99390366 1 495444769693640+07  - . 5938141 90953362 1 29249250+09 
688871331916282943996250+07  - . 108150973332967940126570+09 
594946375302045654710840+07  0. 

341161 465849314922794790+09  . 2706 1 686225238 1 90672825D- 1 5 

598303070674215974075210+08  0. 

0. 


0. 
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.  99390366 149544476969365D+07 
,  68887133 191 628294399627D+07 
.  5949463753020456547 1082D+07 

.  3773 1 994055 1 1 3642365548D+09 


.  845137326302656093460320^-07 
.735713616470395375250780+06 

.  109547062524766292574250+09 


.  845 1 37326302656093460320+07 
.  363955977031654922140340+07 


.  344492538 1 27707978 1 2 1 040+08 
.  151246831629407482626530+06 


.  838931857669547098 1 8 1 720+07 

.971 135293834093658572570+08 
.  230904664605554632432290+08 


.  892952570883608272805980+07 
.  588800794726249031436710+06 
.698325075906120980827600+07 
.  100332698330587598222870+09 
.99301 1151124428630500870+07 


.  58880079472624903 1 436280+06 
.  698325075906 1 209808276 1 0+07 

.239407835148501561946780+09 
.  820010420051903823882810+08 


- . 741907594509337914412860+07 
196332294748688461305000+07 
. 646897926507694997050270+07 


0. 

0. 

- .  679576720894324606539290+08 
- . 168700730230814909976910+09 
- . 309246097320693078172000+08 
0. 

.555111512312578270211750-15 

0. 

0. 

0. 

0. 

0. 

- . 108078477607361602237750+09 
188279777585003039782080+09 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

- .3150727254889985545711 10+08 

-  .  531814122471940895337840+08 

0. 

0. 

- . 684777494941628047641030+07 
0. 

0. 

0. 

0. 

0. 

0. 

. 288100162 1 73630162946960+07 
.714019219249916393789540+07 
0. 

. 122245020841303506087570+07 
.489073674419738001 179440+06 
0. 

0. 

0. 

0. 

0. 

- . 834706847 1 65487480130600+07 

-  .  828325242561074145659830+07 

0. 

. 4890736744 1 973800 1 1 79050+06 
0. 

0. 

0. 

0. 

0. 

-  .  726337294924352223713790+07 
. 256943303904437162802830+08 

-  .  7827858525 1 872984238 1 360+07 
0. 

.395157620137263264751930+07 
. 147063131198997811254030+07 
0. 

0. 

0. 

0. 
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.  90082825256647662065882D+08 
.  18216269902586232160486D+08 


,741907594509337914412850+0? 
.  196332294748688461305000+07 
.  646897926507694997050280+07 

.601690773684571841984590+09 
.  195438529538936275438900+09 


.  7929380619692614  19802060+07 
.  600039355500026774020020+07 
.  50734704 1 700658469 1 58620+07 
.  179040568550176442094690+09 
.538484199201501985566690+08 


.  79293806 1 96926 1 4 1 9802050+07 
.600039355500026774020010+07 
.  50734704 1 700658469 1 58620+07 

.  1 10965020140284615099970+10 
.  322535630794883195961670+09 


.  93 14369 1 7703856331 34 1 700+07 
.  140741659715497925809660+08 
.  276981279912978685734500+07 
.  32020521  1931576805571470+09 
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931436917703856331341700+07 

140741659715497925809670+08 

276981279912978685734490+07 
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0. 

0. 

0. 

0. 

0. 
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- . 736797082928560287759010+07 
0. 

. 108607713702833740073580+08 
. 369062086806381768089810+07 
0. 

0. 

0. 

0. 

- . 224465884074421026229160+08 
•. 175345042843452361656410+09 
- .627257037543626699961560+08 
0. 

. 36906208680638 1 768089810+07 
0. 

0. 

O, 

0. 

0. 

- . 924785042271716594823970+07 
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0. 

.220823718448070651887990+08 
. 667370530455472386000860+07 
0. 

0. 

0. 

0. 
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130391856024987545825480+09 

0. 

. 667370530455472386000880+07 
0. 

0. 

0. 

0. 

0. 

-. 195731365809382768172540+08 
- .  1  17594485280322669256440+09 
- . 368560033447098460246910+08 
0. 

. 353517694054280536201560+08 
. 9883629368 1 8232973702420+07 
0. 

0. 

0. 

0. 
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- . 1 1600235407690759332388D+08 
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0. 
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0. 

0. 

0. 
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0. 
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0. 

0. 

0. 

0. 

0. 

- . 610936065426490388044340+08 
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0. 
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0. 
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0. 

0. 

0. 
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-. 187144983693230792087740+09 
0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

- .317492814674387902095460+08 
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0. 

0. 

- . 583446440083283272440020+08 
0. 

0. 

0. 

0. 

0. 

0. 

. 225489825150449335117390+08 
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0. 
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0. 

0. 
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. 39 16567052 1 3869642 1 68 1 7D+07 

0, 

. 1 19455809396759545039340+10 
. 360868379682357546293790+09 


.  104566988278757288831770+08 
,  243424254361489555579470+08 
.  864279030739663702250440+06 
.  34709405378559 1837308900+09 
.  103008225866985128382570+09 


.  104566988278757288831770+08 
.243424254361489555579480+03 
.  864279030739663702250440+06 

.  104590332773934144733930+10 
.  205502998928787928975730+09 


.  1519744148191 1 1992331840+08 
.255436124896865539205600+08 
.  188884046733634419529190+07 
.304141738363176763799120+09 
.576216669658658653191950+08 


.  15197441481911 1992331840+08 
.  255436 1 24896865539205590+08 
,  188884046733634419529190+07 
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124751087124327178573450+08 
825370888 1 77545398444350+07 

10888679839 1 2663337 1 8 1 70+09 


.  124751087124327178573450+08 
.  1  1 1575550356216078305960+08 


.214524615211724451608720+09 
.  455024126448895550299880+05 


“ . 172371222090186568343990+08 
0. 

. 492556942885473100488210+08 
. 151907712334310389614280+08 
0. 

0. 

0. 

0. 

-.621721 420858979600614920+08 
- .458565086835880446957780+09 
- .  151941221364747724606250+09 
0. 

. 151907712334310389614280+08 
0. 

0. 

0. 

0. 

0. 

- . 200764675997927728868 1 80+08 
-.  131978295878362424042270+09 
- . 429457503375764401035100+08 
0. 

.446475224148703721503420+08 
. 892840097 1 24524708674840+07 
0. 

0. 

0. 

0. 

-. 151464234642419267543080+09 
-.4812171 52790390289253250+09 
- . 999693527745789266914320+08 
0. 

. 892840097 1 24524708674830+07 
0. 

0. 

0. 

0. 

0. 

- . 455658347000786044696690+08 
-. 138800814043570423207540+09 
- . 280203830350265071023240+08 
0. 

. 160470887644577071 164860+08 
0. 

0. 

0. 

0. 

0. 

- , 100073571270047320787970+09 
- . 182705737483255883323280+09 
0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

- . 308098866993286548999930+08 
- , 524740 1 5648 1 686264372580+08 
0. 

0. 

-.21 2459004604689727029450+09 
0- 
0. 
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- . 825706986 1 1 389673305929D+08 

0. 

. 2804164254648 1003059550D+09 
- . 585055924447539380604 1 2D+08 


.  82592479684535876283344D+08 
.719199383493481605185700+07 

.  325954879926238995584200+06 
.  144976495844648261 145990+06 


- .8354487989996840071 19670+08 
- . 832694220076461990227780+05 
. 7 1 68202039757 106486 1 8 1 50+07 
.964434546512351937159160+09 
- .213338224612481 106920130+09 


. 105181447242478195423370+09 
•  .  2374  1  59349 18920588007080+08 
). 

.  144976495844648261 145990+06 


“ . 83269422007646 1 990224860+05 
. 7 168202039757 1 06486 1 8 1 60+07 
0. 

. 100318551728006978086450+09 
- . 88895028 18911 96404390880+07 


.  721  185854674642072365640+07 
.  373885998866773218766220+09 
.816166 112776542947198950+08 

.  117708282581101416552860+07 
.477271818217090527141660+06 


.  723567808845541239514420+07 
.  295490335 1 495673526894 1 0+06 
.  708536536932404003690560+07 
.  281508260526827356758230+09 
.  586 1 52232  1 78333890695280+08 


.  247  162 1 70568 1 37533 1 88760+08 
.  2548975904 1 7096956492430+08 
.816695356763405219475770+07 

.477271818217090527141720+06 


- . 72356780884554 12395 1 4420+07 
- . 295490335149567352689350+06 
. 708536536932404003690570+07 
0. 

.975750134563795764424720+08 
. 2369265477846 1 16187986 1 0+08 


-  . 850373878695554006906 1 30+08 
. 105386557730148053492410+09 
•  . 22503 1 59 1 35496 1 330773850+08 
). 

.408591654522151631656260+07 
.  158372218933617603003270+07 


.731861 165562037873217240+07 
.  104143356512056612637210+07 
.679507418832251508212900+07 
.  100506040512858429913310+09 
.9276121937514231 17817860+07 


.914319641822557122537720+07 
.  864689805228603 147275320+07 
.816413267216772565764940+07 

.  158372218933617603003260+07 


-.731861165562037873217240+07 
- . 1041433565 1 20566 1 2637190+07 
. 6795074 1883225 1 5082 1 2920+07 


•  . 2592444715444  16559840250+08 
.  253534976304450143547070+08 

•  . 6990600744439 15347932620+07 


OOOOOO  IIOOOO  IIOOOOO  O  IIOOOO  IIOOOOO  O  1*0000  1*00000 


120 


.  237208 1 3222937409780545D+09 
.  822008346579402324539030+08 


.  760924333847 1 1695648538D+07 
.  34604056708 1 23 1 57445 1 360+07 
.  587422640505808423680320+07 
.901344574904185068537120+08 
.  187276079369202876895910+08 


.760924333847116956485400+07 
.  34604056708 1 23 1 57445 1 370+07 
.  587422640505808423680330+07 

.663807418278554240268110+09 
.  231  163078961476358858390+09 


.853122655846273581707880+07 
.  122470912739852112070820+08 
.  237258458358110153062280+07 
.  199456458994517219927360+09 
.  655289909206403532404030+08 


.  853 1 2265584627358 1 707870+07 
.  1224709127398521 12070820+08 
.237258458358110153062270+07 

.  80964  1702155993265679050+09 
.  182853726017422380150700+09 


.  120369472887559776244390+08 
.  197120555330874555025000+08 
.217478641425830188428620+07 
.  2386 1 857 1 82 1 167029544580+09 
.  5 1 29086877439686 1 8069200+03 


.  120369472887559776244400+08 
.  197120555330874555025000+08 
.217478641425830188428620+07 

.  35 18269 1 820746342443 1 520+09 


.1261771 9948 1 1 99545657520+08 
. 46 1 607 1 53396740565641020+07 
0. 

0. 

0. 

0. 

- .914602026202446921 194080+07 
- . 549140463828742735737730+08 
- . 2210328021996378 1 3978560+08 
0. 

. 461607153396740565641040+07 
0. 

0. 

0. 

0. 

0. 

-.1041 35256793468525769490+08 
- .792719922766450597698130+07 
- .676224057822054384210290+07 
0. 

. 365349097578288937 1 10550+08 
. 128780518229920044955210+08 
0. 

0. 

0. 

0. 

- . 231039993394594832860020+08 
- .202394459620197855796180+09 
- .797508013146920364208030+08 
0. 

. 128780518229920044955210+08 
0. 

0. 

0. 

0. 

0. 

- . 101906551557583334233570+08 
- . 565550430937343006632440+08 
- .2218610681512554142441 10+08 
0. 

. 452336262720195276570050+08 
. 103455385454325345702820+08 
0. 

0. 

0. 

0. 

- . 808 19167121 296385798 1 980+08 
- . 326594443301215718636670+09 
- . 830 1 295 12871 945656058 1 50+08 
0. 

. 1O34553854543253457C2820+O8 
0. 

0. 

0. 

0. 

0. 

- . 256342304867362068697330+08 
-.941 316785999052572780500+08 
- . 231160713490872622330820+08 
0. 

. 198417311 101749724162160+08 
0. 

0. 

0. 

0. 
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*. 12233586263309947995487D+08 
'  . 868496882401 22953557069D+07 
). 

. 1042C053004015336334178D+09 


12233586263309947995487D+08 

115888149778584492018590+08 


. 704655444503252389679660+09 
• . 689092099 102328928468540+04 


.  906873 153213615921 448730+09 
.  20024 14759 1 2699853 1 30490+09 


.857246279911251727431190+04 


.317875959958337191980040+10 
.  7054887 1 077722 1 297475890-^09 


•  . 2686574017 13319238504940+09 


.  264376598064337391913170+09 
.  544591953394226901634640+08 


.  722002945 1 1067095S46065D+07 
,  7 1 7005030584333954263990+07 


. 90763549 1 372459000960840+09 
199337941933365223206090+09 


.345327522676170224392470+09 


.957097533044302993691670+08 
■ . 719803336296502950067880+07 


- . 84062256778 1 758 1 59460610+08 
167943297188875178621450+09 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

- , 265592151808492755732930+08 
- . 484377849505058160358260+08 
0. 

0. 

- . 703878747923145316226340+09 
0. 

0. 

0. 

0. 

0. 

0. 

,719218598569878824426700+07 

0. 

0. 

. 142976883996446549524140+06 
, 10071887566721 1067467670+05 
0. 

0. 

0. 

0. 

0. 

. 71838983543738476 1 545200+07 
0. 

0. 

. 100718875667211067466280+05 
0. 

0. 

0. 

0. 

0. 

. 121366273386463229120410+10 
0. 

0. 

0. 

. 4008 1 945328405579 1 410900+06 
. 102662558386671832791750+06 
0. 

0. 

0. 

0. 

176748029612390803024460+05 

0, 

0. 

0. 

. 102662558386671822791700+06 
0. 

0. 

0. 

0. 

- . 27034045672057 1 5884 1 5080+09 
- . 7564 1 4039 1552321 90447480+08 
0. 

0, 

0. 
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1 3083300004089 128430887D+06 


.263008517225172794571000+09 
.  5378573364 144075 1 8000820+08 


.773267656348750520801810+08 
.211 795375362238000852640+08 


.  96677708709075455 1731 280+08 
.  245740521473753009432490+08 


.728010299255097061169210+07 
.  69723205572526575988 1120+07 


.  94970106368201 1946846580+08 
.  80800 1 506022272206 1 1 2440+07 


.  22330972634 1 985996 1  1 5630+08 


.  283440375017551667661550+09 
.  1121516727908867771  18280+09 


- .2182229494401 16282106620+07 


.  101037250125997786560650+09 
.  289426920316660874020020+08 


.  9 1 76758050744092058 1 1 100+07 
.  997998603770865614766390+07 


.497472446136152583855990+09 
.  140375599292062531563270+09 


.  82 1  16093765 1 726240123700+07 
.  54804 14799 1 4675079 1 767 10+07 


- . 7234 1 1 376999087 1 24584300+07 
.712490939216188327583250+07 
0. 

0. 

0. 

. 4373706448 1 8506297 187930+06 
0. 

0. 

0. 

0. 

.961 84648905 1 3000390024 1 0+08 
0. 

0. 

0. 

0. 

. 4075 1 98903 1 3593560624580+07 
. 151069223487681164712640+07 
0. 

0. 

0. 

-.518921 30658040405821 5630+06 

0. 

0. 

0. 

o. 

. 151069223487681164712620+07 
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